ABSTRACT 


A  mathematical  model  of  a  pyrotechnic  flame  has  been 
developed  and  validated.  The  model  utilizes  a  computer  program 
to  calculate  the  relative  radiant  pov/cr  of  flares  containing  any 
of  the  alkali  metals.  The  computation  takes  into  consideration 
known  system  variables  such  as  formulation,  candle  size, 
displacement  along  the  plume  axis,  and  atmospheric  pressure. 

These  effects  were  evaluated  over  a  wide  range  of  conditions  and 
the  computed  spectra  follow  the  established  experimental  spectra 
trends.  The  computer  solution  is  a  general  solution  requiring  no 
ad  hoc  modifications  to  -'the  model  for  any  of  the  metals  or 
conditions  tested. 

The  basis  for  the  model  of  the  flame  permits  calculation  of 
the  effects  of  mixing  of  the  flare  plume  with  the  ambient 
atmosphere  and  takes  heat  loss  along  the  plume  axis  into 
consideration.  The  relative  radiant  power  is  determined  by 
solving  the  radiative  transfer  equation  using  parameters 
calculated  from  the  thermodynamic  properties  of  the  flare  plume. 
The  individual  line  dispersion  profiles  of  the  resonance  lines 
were  replaced  with  a  two-line  Voigt  function  that  includes  a 
broadening  parameter  that  considers  natural,  Lorentz,  resonance, 
quenching  and  Doppler  broadening  mechanisms.  The  solution 
results  in  a  spectral  profile  that  follows  observed  experimental 
spectral  variations  over  a  wide  range  of  metals  and  conditions. 
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Spectra  were  computed  for  lithium,  sodium,  potassium, 
rubidium  and  cesium  under  a  variety  of  limiting  conditions, 
comparison  of  measured  spectra  v/i  th  the  computed  spectra  showed 
good  agreement  both  qualitatively  and  quantitatively.  Cases 
shown  include  theoretical  and  experimental  spectra  comparing 
varying  ambient  pressure,  formulations,  displacement  along  the 
plume  axis  and  candle  diameters  for  sodium  and  selected  cases  for 
the  other  alkali  metals  to  compare  with  experimental  spectra  or 
show  expected  trends. 

This  model  should  be  useful  in  optimizing  flare  radiant 
output  using  any  of  the^alkali  metals  as  emitters  and  also  as  a 
tool  for  further  study  of  the  radiative  transfer  system  in 
controlled  laboratory  conditions. 
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INTRODUCTION 


Early  research  in  magnesium-sodium  nitrate  flame 
spectroscopy  has  shown  that  the  major  component  of  the  visible 
light  emission  from  these  flames  is  a  continuous  spectrum  called 
the  sodium  resonance-line  continuum.  The  resonance-line 
continuum  is  much  broader  than  predicted  by  Doppler  or  Lorentz 
effects.  Further,  many  observations  of  flame  spectra  seem  to 
differ  in  the  effects  noted  due  to  varying  test  parameters. 

These  incongruities  led  researchers  such  as  Douda,  Blunt  and  Bair 
(33,48,49)  to  consider  other  mechanisms,  such  as  radiative 
transfer,  to  elucidate  the  origin  of  the  resonance-line 
continuum.  Their  early  work- indicated  that  radiative  transfer 
was  able  to  yield  reasonable  agreement  in  specific  cases  but 
required  certain  ad  hoc  modification  to  fit  various  flame  cases. 
Douda  (1)  hypothesized  that  radiative  transfer  was  the  major 
mechanism  responsible  for  this  resonance-line  continuum.  He 
developed  equations  for  predicting  broadening  of  the  sodium  D- 
lines  using  some  simplifying  assumptions  for  temperature  and 
concentration  profiles  and  neglecting  the  effect  of  air  mixing 
with  the  flame.  By  comparing  the  spectra  of  flames  at  reduced 
pressure,  he  demonstrated  that  the  theory  produced  reasonable 
agreement  with  the  experiments.  The  agreement  was  improved  at 
the  lower  pressures  because  the  effect  of  the  atmosphere  mixing 
with  the  flame  was  reduced.  Even  with  this,  however,  the 
coefficient  used  in  the  spectral  line  profile  had  to  be  reduced 
by  almost  an  order  of  magnitude  to  get  agreement  between 
calculations  and  experiments.  In  1975  and  1976,  the  model  was 
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expanded  by  Dillehay  (4,5)  to  take  into  account  the  effects  of 
air  nixing  with  the  flame.  This  work  improved  the  calculation  of 
the  temperature  profile  to  obtain  agreement  with  experimental 
observations.  The  sodium  atom  concentration  was  not  assumed  to  be 
a  constant  average  value  in  the  flame  but  was  approximated  as  a 
parabola  through  the  flame  cross-section.  Using  these  data,  the 
line  broadening  profiles  were  calculated  for  a  number  of 
conditions  without  making  any  a  priori  assumptions.  The 
predicted  spectra  showed  the  same  broadening  fe."  res  as  the 
experimental  data  even  at  atmospheric  pressure.  'hese 
calculations  were  all  done  by  hand  with  considc  le  use  of 
approx ima ti ons ,  interpolations  and  look-up  tabli  ■  The  difficult 
and  time  consuming  hand  calculations  in  this  early  model  limited 
the  application  of  this  technique  for  further  research. 

The  purpose  of  this  dissertation  will  be  to  develop  an 
improved  mathematical  flame  model  that  can  be  used  to  verify  the 
radiative  transfer  theory  for  sodium  and  add  other  alkali  metal 
species  to  the  model.  By  generating  a  general  solution  to  the 
equations,  the  model  will  allow  computation  of  the  broadened  line 
profiles  of  lithium,  sodium,  potassium,  rubidium  and  cesium.  The 
solution  will  be  a  computer  simulation  of  a  pyrotechnic  flame  and 
a  computer  program  that  will  take  the  predicted  combustion 
species  from  the  simulation,  calculate  the  radiative  transfer 
parameters  and  generate  the  spectral  profile.  The  flame  model 
will  take  into  account  the  effects  of  varying  formulation,  mixing 
of  the  atmosphere  with  the  flame,  heat  loss  along  the  flame  axis 
and  the  effects  of  pressure  on  the  flame. 


The  equilibrium 
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compositions  and  temperatures  are  calculated  using  the  NASA  SP- 
273  Combustion  Equilibrium  Program  (3).  This  model  will  enhance 
the  state-of-the-art  of  illuminant  formulating  and  improve  the 
understanding  of  the  radiative  energy  transfers  in  pyrotechnic 
flames. 

This  research  will  show  that  the  radiative  transfer  model 
can  be  applied  to  pyrotechnic  flames  without  a  priori  assumptions 
or  ad  hoc  modification  of  the  model.  The  model  will  then  be 
used  to  demonstrate  that  many  of  the  variables  noted  in  previous 
experiments  are  predicted  by  the  theory  and  quantifying  rules  can 
be  established  for  the  meaningful  acquisition  and  interpretation 
of  pyrotechnic  flame  spectra.  The  spectral  model  will  take  into 
account  the  effects  of  varying  formulation,  candle  diameter, 
pressure  and  position  within  the  flare  plume.  Data  for  lithium, 
sodium,  potassium,  rubidium  and  ces-uin  will  be  included  to  allow 
a  broader  application  of  the  model  and  increase  the  confidence  in 
the  validation.  Experimental  data  are  already  available  for  a 
variety  flare  formulations  taken  under  many  conditions.  These 
data  will  be  used  to  compare  with  the  calculated  spectral 
profiles. 
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THEORETICAL  BACKGROUND 


The  Radiative  Transfer  Equation 


A  detailed  derivation  and  discussion  of  the  radiative 
•transfer  equation  can  be  obtained  from  the  literature  (1,43,50). 
A  brief  summary  will  be  presented  here  to  aid  in  the 
understanding  of  the  model. 


To  characterize  the  mechanisms  of  luminous  energy  generation 
in  a  flare  plume,  it  is  necessary  to  solve  the  general  radiative 
transfer  equation. 


\  =  Exp [-(r-r, ) 


r-r, 

where  TJ  and  rz  are  the  optical  depth  integration  limits  from 
front  to  rear  of  the  plume  respectively ,  T^j  and  are  the 

spectral  intensities  at  optical  depths  ZJ"  and  7^  ,  S^Ct)  is  the 
line  source  function,/?  is  the  cosine  of  the  angle  of 
observation,  and  YveL  is  the  normalized  spectral  profile  of  the 
absorption  coefficient.  To  solve  this  equation,  the  plume  is 
constrained  by  the  following  boundaries: 


Ci) 


1)  The  plume  is  a  homogeneous  gaseous  atmosphere  with  plane- 
parallel  stratification. 


2) 


The  gas  consists  of 

excited  to  the  ZRy 

'z 

which  can  interact 


alkali  metal  atoms  which 
or  R,  level  and  molecular 
to  perturb  the  metal  atom 


can  be 
species 
radiation.  . 


3)  There  is  local  thermodynamic  equilibrium  (LTE)  governed 
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by  the  local  temperature. 

4)  Energy  exchange  by  radiation  leads  to  radiative 
equilibrium, 

5)  The  refractive  index  of  the  medium  is  unity. 

6)  The  radiation  is  unpolarized  when  emitted  and  remains 
unpolarized  on  its  interactions  with  flare  species. 

7)  The  temperature  gradient  can  be  constructed  from  the 
thermodynamic  equilibrium  temperature  calculated  for  a 
limited  surface  penetration  of  the  ambient  atmosphere. 

8)  The  absorption  profile  ^  and  the  number  density  of 
sodium  atoms  are  dependent  on  the  calculated 
thermodynamic  equilibrium  compositions  defined  by  the 
limited  surface  interactions  with  the  ambient  atmosphere. 

Using  these  boundary  conditions,  the  radiative  transfer  equation 
can  be  simplified  by  assuming  that  the  observed  flux  is  emerging 
normal  to  the  surface  {(i  =  1)  and  there  is  no  flux  incident  on 
the  rear  surface  of  the  plume  (1^  =  0)  and  for  the  LTE  case 
Sm(t)  -  B„(T')  where  is  the  Plank  function  for  the  flare 

temperature  at  optical  depth  T  .  The  monochromatic  emergent 
intensity  is  now 

t-T 

I„  s  C-r$Jdr 

where  T  is  the  total  optical  thickness. 


(2) 
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A  numerical  integration  of  Eq.(2)  using  Simpson's  Rule  of  2M 


intervals  (2)  yields  a  theoretical  relative  spectral  radiant 
power  for  a  given  wavelength,  proportional  to  spectral  emergent 
intensity.  By  calculating  values  over  a  spectral  frequency 
range,  a  theoretical  flare  spectrum  can  be  generated. 


In  order  to  solve  the  transfer  equation  simultaneously  for 
both  of  the  sodium  D-lines,  Douda,  Behmenberg,  Hoffman  and  Kohn 
(1,21,26,27)  replaced  the  individual  dispersion  profiles,  in 

Eq. (2) ,  with  a  2-line  Voigt  function.  This  function  includes  a 
superposition  of  broadening  mechanisms  and  is  described  as 


<90 


where  V  is 
of  Doppler 


cl  Gyp  (~yz) 
(y-y)z*d-  W 


— o o 

the  frequency  measured  from  the  line  center 
half-width  and  a.  is  the  dispersion  parameter 


in 


(2) 

uni  ts 


The  dispersion  parameter,  a.  ,  is  obtained  from  the  equation 

y 

6L  -  h-aAl*+-  AXf>  +  AA^/AX^]  2.)  (4) 

where  aAn  ,  AAl  ,  AAR  ,  aA^  ,  and  are  the  natural,  Lorentz 

(unlike  particle),  resonance  (like  particle),  quenching,  and 
Doppler  line  broadening  half-widths,  respectively. 


The  Voigt  Function  0.  Parameter 

The  Voigt  <X.  parameter  shown  in  Eq.(3)  is  the  dispersion 
parameter  that  takes  into  account  all  of  the  factors  considered 
that  contribute  to  the  width  of  the  spectral  line.  The  line 
broadening  mechanisms  considered  in  the  pyrotechnic  flame  spectra 
are  (a)  natural  broadening  (due  to  the  Heisenburg  uncertainty 
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principle  (44)),  (b)  Holtsmark  pressure  broadening  due  to 
collisions  with  like  neutral  species,  (c)  van  der  Waals  pressure 
broadening  due  to  collisions  with  unlike  neutral  species,  (d) 
broadening  due  to  quenching  collisions,  and  (e)  Doppler 
broadening  due  to  the  relative  motion  of  the  radiating  system  and 
the  observer.  Stark  broadening  due  to  charged  perturbers  was 
neglected  as  the  concentration  of  charged  metal  ions  in  the  flame 
can  be  shown  to  be  small  (D)  .  As  shown  in  Eq.(4),  the  effect  of 
all  the  considered  broadening  mechanisms  is  accounted  for  in  the 
Voigt  function  a.  parameter.  Each  of  the  line  broadening  half¬ 
widths  will  be  discussed  separately  below. 

The  equation  for  natural  line  broadening  is 

A\(A)  ~  lo8* o  A2,/?7rc  & 

where  A0 is  the  line  center  wavelength  and  is  the  Einstein 
coefficient  for  stimulated  emission.  The  values  of  Azl  for  the 
alkali  metals  are  shown  in  Table  1.  along  with  the  calculated 
natural  line  broadening  for  each. 

Table  1. 

Natural  Line  Broadening  and  Einstein 
Coefficients  for  the  Alkali  Metals 


Metal 

,  -Sfic 

Li 

6707 

. 37E+08 

0.9E-04 

Na 

5893 

. 65E+08 

1.2E-04 

K 

768  2 

.  39E+08 

1. 2E-04 

Rb 

7874 

. 36E+08 

1.2E-04 

V 
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The  Doppler  broadening  half-width  for  the  alkali  metals  is 
given  by 

AX0(A)  =  2xio8[C2<RA2)^/c]^(T0/M)y2  (6) 

where  M  is  the  molecular  weight  of  the  metal  in  gm./mole.  Table 
2.  shows  the  values  calculated  for  the  Doppler  broadening  half¬ 
width  for  a  temperature  of  3000  °K. 

Table  2. 

Doppler  Line  Broadening 
for  the  Alkali  Metals 


Metal 

Xo 

M 

A 

Li 

6707 

6.94 

0.0999 

Na 

5893 

22.99 

0.0482 

K 

7682 

39.10 

0.0482 

Rb 

7874 

85.48 

0.0334 

Cs 

8732 

132.91 

0.0297 

Since  the  value  of  is  proportional  to  the  square  root  of  the 

flame  temperature ,  aXp  for  sodium  varies  from  0.0482  A  at  3000  °K 
to  0.0387  A  at  1915  °K. 

The  van  der  Waals  collisional  broadening  (also  known  as 
Lorentz  broadening  )  due  to  unlike  species  is  given  by 

aXl(A)  =  /oe0\t/c)<rLN*  (8kX,^A/TT^)y2  (l) 

where  <5^,  is  the  optical  cross  section  (cnr  )  of  the  alkali  metal, 
and  is  the  number  density  of  perturbing  species  (cm  ).  is 

Avagadro's  number  and  f*-  =  mM/(m+M),  the  reduced  mass.  If  the 
mean  relative  velocity  of  colliding  species  (cra'Sec'1 )  is 

--  C8kT.NA/-v.>*  ®> 

where  j  is  the  species  index , . Eq.  (7)  can  be  written 

a\(&)  =  io6()£/c)  <rL  ^ 
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The  values  for  optical  cross  section  at  2500  °K  reported  by 
Hoffman  and  Kohn  (21)  are  shown  in  Table  3.  below,  along  with 
typical  values  for  the  Lorentz  broadening. 

Table  3. 

Optical  Cross  Sections  and 
Lorentz  Line  Broadening  Values 
for  the  Alkali  Metals 


Metal 

< TuCxio 

A^l.  j  A 

Li 

46.5 

0.05880 

Na 

59.3 

0.02967 

K 

59.1 

0.04043 

Rb 

76.3 

0.03313 

Cs 

89.1 

0.04355 

The  values  for  Lorentz  broadening  half-width  in  Table  3.  do  not 
take  into  account  the  effects  of  mixing  the  ambient  atmosphere 
with  the  plume.  These  effects  will  be  considered  at  the  end  of 
this  section. 

The  Holtsmark  or  resonance  broadening  half-width  can  be 
derived  from 

A^CA)  =  [C6*loVM  X?  ]  /r^7T2c2  me)  (10) 

where  e  is  the  elementary  charge  (cm  -  gm  sec  )  ,  is  the 
electron  rest  mass  (gm)  ,and  A/j  =  is  the  number  density  of  the 
alkali  metal  atoms.  The  oscillator  strengths,^,  of  the  alkali 
metals  were  obtained  from  several  sources  (39,40,42).  The 
uncertainty  of  the  determined  values  is  not  of  sufficient 
magnitude  to  cause  any  significant  variation  in  the  calculations. 
Resonance  broadening  of  the  spectral  lines  is  the  largest  factor 
in  determining  the  Voigt  function  Q.  parameter  except  at  very  low 
metal  atom  concentrations.  Again  without  considering  the  effects 
of  interactions  with  the  ambient  atmosphere,  the  values  of 
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resonance  broadening  half-widths  and  oscillator  strengths  are 
shown  in  Table  4. 


Table  4. 

Oscillator  Strengths  and 
Resonance  Line  Broadening  Values 
of  the  Alkali  Metals 


Metal 

f, 

ft 

/ 

;  A 

Li 

0.502 

0.251 

0.753 

0.52342 

Na 

0.624 

0.312 

0.936 

0.63442 

K 

0.682 

0.339 

1.021 

1.30661 

Rb 

0.675 

0.322 

0.997 

1.19584 

Cs 

0.732 

0.362 

1.09 

1.88006 

The  equation  for  quenching  line  half-width  can  be  expressed 


as 

aXq(?\)  *  /o8A*  C2iCt)/2ttc  (ll) 


where  is  the  rate  of  collisional  deactivation  (  i.e.,  the 

number  of  quenching  collisions  per  second).  Hooymayers  (19,20) 
equates  the  quenching  process  to  a  shortening  of  the  radiative 
lifetime.  He  calculates  the  value  of  CZl(T)  as 

c2i(t)  ^  qz1ct)  *  nj  vj  Qz) 

where  vij  is  the  density  of  the  quenching  flame  specie  j ,  Vj  is 
the  average  relative  velocity  of  approach  of  the  colliding  specie 
j,  and  Gj  is  the  specific  quenching  cross  section.  When  the 
specific  alkali  metal  atoms  undergo  quenching  collisions  in  a 
mixture  of  flame  molecules  of  different  species,  the  quenching 
frequency  is  given  by 


where  j  is  the  index 
The  value  for  Vj  is  g 

*4 


for  the  various 
iven  by 


species 

K 


of 


flame  molecules. 


Page  11 

The  mean  speed  of  the  alkali  metal  atom  and  of  any  specie  atom  j 
is  given  by 

'W’V.  «  CSkX^ArM,)-*  (is) 

where  is  Avagadro's  number  and  Mj  is  the  molecular  weight  of 
the  metal  or  specie  j.  The  ideal  gas  law  is  used  to  calculate 
the  number  density  of  the  species  in  the  flame  by  relating  the 
mole  fraction  of  the  specie  j  to  the  total  gas  fraction  at  a 
given  temperature  and  pressure.  Values  for  the  quenching  cross 
section  for  various  species  with  the  different  alkali  metals  were 
required  for  the  computation  of  Cit£T)  •  Hooymayers  and  Alkemade 
(19,20)  reported  measured  values  for  quenching  cross  sections  of 
various  species  with  sodium  and  potassium.  It  is  thought  that 
the  apparent  quenching  cross  section  of  a  specie  with  an  alkali 
metal  would  be  inversely  proportional  to  the  atomic  weight  of  the 
metal.  This  would  be  the  result  of  the  collisional  mass  of  the 
metal  varying  while  the  mass  of  the  quenching  specie  remained 
constant.  The  quenching  cross  sections  of  species  in  collision 
with  sodium  have  been  examined  more  extensively  than  the  other 
alkali  metals.  The  values  of  quenching  cross  sections  with  sodium 
were  divided  by  the  ratio  of  atomic  weights  of  potassium  and 
sodium  (39/23)  to  predict  the  values  for  potassium.  The 
calculated  values  were  compared  with  the  measured  values  of 
Hooymayers.  The  results  are  shown  in  Table  5. 


I 
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Table  5. 


Quenching  Cross  Section  Values  for 
Species  with  Sodium  and  Potassium 


0j  Cx  lO>0) /c^2- 

Specie 

Temp 

N a  Value 

Calc.  I<  Value 

Meas . 

K  Value 

CO 

2000 

72 

42 

40 

C02 

1670 

113 

66 

80 

Ar 

2070 

2.3 

1.35 

0.94 

H20 

2070 

1.0 

0.59 

0.42 

N2 

1940 

40 

23 

18 

02 

1885 

66 

39 

31.5 

The  agreement  is 

entirely  satisfactory  given 

the  unci 

ertainty  o 

f 

the  original  mea: 

surements.  Indeed,  variations 

g  reate 

r  than  the 

se 

are  noted 

betweei 

n  results  from 

different  rese 

archers 

.  The 

quenching 

cross 

section  values 

of  all  of  the 

major  s 

pecies  in 

the 

flare  flame  were 

calculated  with  reference  to 

the  sodium 

quenching 

cross 

section  values 

and  the  atomic 

we  i  g  h  t 

of  the 

alkali  metal.  These  values  are  tabulated  in  Table  6.  below  for 
2  0  0  0  °K . 


Table  6. 

Quenching  Cross  Section  Values  of  Species 


with 

the  Alkali  Metals 

C*I0I<5),C 

Specie 

Li 

Na 

K  Rb 

Cs 

Ar 

7.9 

2.4 

1.4  .6 

.4 

CO 

239 

72 

42  19 

12 

C02 

312 

94 

55  25 

16 

H 

3.3 

1 

.6  .3 

.2 

112 

53 

16 

9.4  4.3 

2.8 

1120 

14.6 

4.4 

2.6  1.2 

.8 

Mg 

3.3 

1 

.6  .3 

.2 

N  2 

139 

42 

25  11 

7.3 

02 

205 

62 

36  17 

11 

These  quenching  cross  sections  were  used  in  the  computation  of 
the  quenching  broadening  half-width.  Typical  values  of  the 
quenching  broadening  half-width  calculated  for  each  metal  without 
regard  for  the  interaction  with  the  air  are  shown  in  Table  7. 
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Table  7. 

Quenching  Lino  Broadening 
of  the  Alkali  Metals 


Metal 

t  A 

Li 

4 . 69  54E-0  3 

Na 

1. 9329E-03 

K 

1 . G166E-0  3 

Hb 

1.0645E-03 

Cs 

6.2197E-04 

From  Eq.(4),  the  Voigt  function  cl  parameter  can  be 
calculated  for  each  metal.  Without  taking  into  consideration  the 
effect  of  the  ambient  atmosphere,  the  values  for  the  ex  parameter 
for  stoichiometric  formulations  are  shown  is  Table  8. 


'  Table  8. 

Voigt  Function  CL  Parameter  Values 
for  the  Alkali  Metals 


tal 

a~ 

Li 

4.8372 

Na 

11. 3G22 

K 

23.1741 

Rb 

38.8102 

Cs 

54.7553 

To  examine  the  effect  of  mixing  the  ambient  atmosphere  with 
the  flare  plume,  thermochemical  calculations  were  made  with  the 
admixture  of  air  with  the  flare  composition.  Calculations  were 
run  in  10%  increments  from  0%  air  to  90%  air.  From  these  data, 
new  values  of  the  Voigt  function  parameter  were  calculated. 

Table  9.  and  Figure  1.  show  the  values  obtained  for  sodium  over 
the  range  of  formulations. 


•  t 
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Table  9. 

Line  Broadening  and  Voigt  Function  ct  Parameter  Values 


for 

Sodi urn 

as  a  Function  of 

Z  Air 

Admi xed 

%Ai  r 

T,°k 

A  X  u 

aXl 

AX* 

aXq 

CL 

0 

3078 

.00012 

.02967 

.  63442 

.00193 

.04881 

11.362 

10 

3064 

.00012 

.03144 

. 54871 

.00216 

.04070 

9.9  57 

20 

3038 

.00012 

.03303 

.46972 

. 00238 

.04849 

8.674 

30 

2999 

.00012 

.03453 

. 39649 

.00257 

.04818 

7.495 

40 

2938 

.00012 

.03556 

. 32874 

.00274 

.04769 

6.410 

50 

2835 

.00012 

.03632 

. 26576 

.00280 

.04684 

5.422 

60 

2630 

.00012 

.03632 

. 20560 

.00295 

.04512 

4.521 

70 

2246 

.00012 

.03491 

.13779 

.00287 

.04170 

3.508 

80 

1710 

.00012 

.03135 

. 07827 

. 00260 

.03638 

2.571 

90 

1181 

.00012 

.02662 

. 00006 

.00199 

.03023 

0.793 

From  the  table,  it  can  be  seen  that  the  Doppler  half-width  is 
only  reflecting  changes  in  the  temperature  of  the  mixtures.  The 
natural  broadening  half-width  is  only  a  function  of  the 
wavelength  and  the  Einstein  coefficient  of  stimulated  emission. 
The  Lorentz,  resonance  and  quenching  parameters  are  influenced 
not  only  by  the  temperature  but  also  by  the  reduced  sodium  atom 
concentration  and  the  change  in  the  perturbing  species.  It  is 
particularly  notable  that  the  nitrogen  and  oxygen  content  of  the 
plume  mixture  increase  markedly  since  they  are  especially 
effective  in  quenching  sodium  D-line  radiation.  All  of  the 
alkali  metal  Voigt  function  a  parameters  show  a  similar  response 
to  the  mixing  of  air  with  the  flame.  Again  using  sodium  as  an 
example,  Table  10.  and  Figure  2.  show  the  effect  of  varying  the 
formulation  of  the  illuminant  on  the  Voigt  function  Q.  parameter 
at  the  0%  air  level. 
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Table  10. 

Line  Broadening  and  Voigt  Function  Parameter  Values 


for  S< 

ad i  urn  as  a 

Function 

of  %  ! 

Mngnesi urn 

at  81 

Bi  ndc  r 

%Mg 

T/k 

A  Xf> 

<X 

30 

294  0 

.00012 

.02425 

.79960 

.00159 

.04801 

14.317 

35 

2900 

.  00012 

.02543 

.73523 

.00144 

.04845 

13.098 

40 

2953 

. 00012 

.02761 

. 65122 

.00133 

.04824 

1 1 .742 

4  5 

2015 

. 00012 

.03016 

. 56543 

.00126 

.04759 

10.443 

50 

2350 

.00012 

.03315 

.  49993 

.00123 

.04606 

9.660 

55 

1964 

.00012 

.03817 

. 48608 

.00125 

.04201 

10.417 

60 

1938 

.00012 

.04251 

.48479 

.00106 

.03893 

11.303 

65 

1795 

.00012 

.04275 

.  47700 

.00051 

.03830 

11.313 

It  may  be  seen  that  the  Voigt  function  CL  parameter  is  sensitive 
to  the  changes  in  formulation.  While  only  sodium  values  are 
shown  in  the  table,  all  of  the  alkali  metals  show  a  similar 
response  to  changing  formulation. 


Optical  Thickness 


To  evaluate  Eq.(2),  a  means  is  needed  to  determine  the  total 
optical  thickness  for  a  given  point  within  the  plume.  Douda  used 
an  equation  for  optical  thickness  that  assumed  that  the  number 
density  of  sodium  atoms  had  an  average  value  inside  the  flame. 
With  this  assumption,  the  optical  thickness  can  be  described  as 

X(?)  =  [\uoATr^yJ(N ,  6(2  -  Uz6zi)Z 

where  and  BZ\  c'iro  the  Einstein  absorption  and  induced 

emission  coefficients,  is  the  number  density  of  metal  atoms  in 
the  lower  { )  state  and  is  the  number  density  of  metal  atoms 
in  the  upper  (ZP  )  state  and  Z  is  the  flame  depth  (cm.).  The 
Doppler  hal f -width ,  Ai/e  ,  is  in  units  of  frequency.  It  can  be 
shown  that  N z  is  negligible  compared  to  (28,41).  Eq.(16)  is 

not  adequate  for  the  model  described  in  this  research  since  the 


metal  atom  concentration  is  not  assumed  constant  but  allowed  to 
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vary  as  a  function  of  position  within  the  plume.  Accordingly, 
the  model  constructs  a  profile  of  the  metal  atom  concentration 
across  a  section  of  the  plume  and  calculates  the  optical 
thickness  by  numerical  integration  using  the  trapezoidal  rule  to 
give  the  value  at  any  depth  in  the  plume.  See  Appendix  B  for  a 
detailed  description.  The  expression  for  optical  thickness  may 
be  written  as 

r  =  rU)  f  07) 

where  is  the  metal  atom  concentration  in  grnm-atoms/cc  at 

depth  Z.  Figure  3.  shows  the  optical  thickness  as  a  function  of 
depth  in  the  flame  for  a  typical  sodium  composition. 
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THE  FLARE  FLAME  MODEL 

Physic;:!  end  Optical  Boundaries 

To  model  the  combustion  of  a  flare  composition  one  must  take 
into  account  several  factors  that  are  known  to  influence  the 
flame  both  physically  and  chemically.  As  the  fuel  and  oxidizer 
react,  the  combustion  products  are  expelled  from  the  surface  and 
form  a  gaseous  volume  that  radiates  both  line  spectra  and  a 
continuum  due  in  part  to  radiative  transfer  broadening  and  in 
part  to  gray  body  emission  from  the  solid  particles  in  the  plume. 
This  research  will  consider  only  the  radiative  transfer  aspects 
of  the  resonance-1 '■  ne  continuum. 

For  this  model,  the  burning  rate  of  the  composition  is 
assumed  to  be  optimized  for  the  formulation.  This  is  attainable 
in  practice  by  tailoring  the  particle  size  distribution  of  the 
magnesium  or  other  fuels.  It  has  been  observed  that  there  is  an 
optimum  burning  rate  for  a  given  formulation.  If  the  rate  is 
slower  than  the  optimum,  heat  losses  will  be  excessive  and  a  low 
light  integral  will  be  recorded.  If  the  rate  is  faster  than  the 
optimum,  material  will  pass  through  the  reaction  zone  too  fast  to 
allow  equilibrium  combustion  and  material  will  be  ejected  from 
the  plume  without  contributing  to  the  system  resulting  in  a  low 
light  integral.  This  model  assumes  that  proper  adjustment  of  the 
burning  rate  has  been  accomplished. 
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The  stoichiometric  formulation  of  a  composition  is  the 
formulation  that  gives  the  highest  temperature  when  the 
equilibrium  proceeds  to  minimize  the  Gibb's  free  energy  of  the 
system.  This  is  usually  a  balance  point  where  just  sufficient 
oxidizer  is  present  to  oxidize  the  fuel  present.  If  the  reaction 
products  are  allowed  to  mix  with  the  ambient  atmosphere,  however, 
and  the  ambient  atmosphere  is  air,  there  is  a  potential  to 
further  react  excess  fuel  with  the  oxygen  in  the  air.  This  will 
proceed  to  a  point  where  there  is  not  sufficient  fuel  reacting 
with  the  air  to  raise  the  temperature  of  the  air  to  the 
combustion  temperature.  As  the  flare  formulation  is  adjusted  to 
be  fuel  rich,  the  excess  fuel  reacts  with  the  air  and  the 
reaction  proceeds  further  along  the  plume  axis  to  give  a  longer 
plume.  This  has  been  observed  experimentally  by  many 
investigators  (11,12,13,14,15,16,24,29,30,34,35,36).  The  plume 
length  is  also  a  function  of  the  pressure  surrounding  the  plume. 
The  reduced  pressure  allows  the  plume  to  expand  along  the  plume 
axis  resulting  in  longer  plumes.  This  is  reported  by  Douda  and 
Blunt  (11,15,16).  The  length  of  the  flare  plume  is  also  a 
function  of  the  diameter  of  the  flare.  As  the  cross  section  of 
the  plume  increases  in  diameter,  it  obviously  takes  longer  for 
the  plume  to  mix  with  the  air  and  complete  the  combustion. 

To  solve  the  radiative  transfer  problem  requires  a 
mathematical  model  of  the  flare  plume  to  permit  computation  of 
the  temperature  and  metal  atom  concentration  profiles  through  the 
plume  at  a  given  height  above  the  flare  surface.  Figure  4.  shows 
the  physical  model  that  was  derived.  The  outer  boundaries  are 
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represented  by  parabolas  that  are  fitted  to  the  edge  of  the 
candle  and  extend  to  a  point  above  the  candle  that  is  defined  as 
a  function  of  the  candle  diameter  and  formulation.  The  parabolas 
are  defined  as  the  902  air  interface  and  represent  the  physical 
boundaries  of  the  plume.  The  factor  F  is  the  fuel  fraction  that 
gives  a  computed  temperature  of  1800  °K.  Lewis  and  von  Elbe  (48) 
reported  that  sodium  D-line  emission  disappears  below  1775  °K. 
Fristrom  and  Westenbcrg  (45)  reported  that  omission  of  additives 
to  flames  ceases  at  about  1700  °K.  The  thermochemical  combustion 
analysis  shows  that  a  limit  of  1800  °K  agrees  well  with  the 
disappearance  of  the  metal  atoms  from  the  plume  due  to  the 
formation  of  the  oxides  and  hydroxides.  Brewer  and  Searcy  (6,7) 
report  that  the  oxides  of  the  various  alkali  metals  are  stable  at 
temperatures  in  the  1500  to  2600  °K  range  with  lithium  oxide 
forming  at  2600  °K  and  cesium  oxide  forming  around  1500  °K.  It  is 
assumed  that  an  average  temperature  of  1800  °K  will  be  suitable 
for  all  the  alkali  metals  for  the  disappearance  of  emission.  The 
factor  F  varies  as  the  formulation  varies  and  for  a  given 
diameter ,  (d)  ,  the  plume  lengthens  as  the  composition  moves 
towards  a  fuel  rich  composition.  Photographs  of  flare  plumes 
(32,36)  show  good  agreement  with  the  values  calculated  for  plume 
length.  The  maximum  height  of  the  visible  plume  above  the  candle 
surface  is  given  as 

-  dL/F  (18) 

To  account  for  the  influence  of  ambient  pressure  on  the  plume 
length,  another  term  is  added  to  the  length  calculation  to  give 

L„„  =  (™o/p)Kj  01) 


Page  20 


The  cube  root  of  the  pressure  relative  to  1  atmosphere  gives 
reasonable  agreement  with  the  observed  spectra  previously  cited. 
The  diameter  (DL)  of  the  plume  at  the  maximum  height  above  the 
flare  surface  is  determined  as 

=  (7i 0/p/6d(U  \/r-)'/z  (go) 

The  adjustment  of  length  and  diameter  for  pressure  allows  the 
volume  of  the  plume  to  vary  with  the  pressure  according  to  the 
Ideal  Gas  Law  relationship  fjvj  =■  &,  \'z  .  These  two  parameters 

permit  the  computation  of  the  outer  boundaries  of  the  visible 
plume.  To  account  for  the  mixing  of  the  ambient  atmosphere  with 
the  plume,  a  parabolic  Inner  boundary  was  established  as  the  01 
boundary  for  admixture  of  the  atmosphere  from  one  direction  into 
the  plume.  The  parabola  used  has  the  same  parameters  for  all 
diameters  of  candles.  This  would  be  expected  as  entrainment  and 
diffusion  should  be  independent  of  the  diameter  of  the  candle. 
This  boundary  is  mirror  imaged  to  account  for  the  admixture  of 
the  atmosphere  from  the  other  side  of  the  plume.  It  was  felt 
that  the  air  would  be  rapidly  introduced  near  the  surface  due  to 
turbulence  and  entrainment  but  that  diffusion  would  be  the 
primary  means  of  transport  within  the  inner  regions  of  the  plume. 
Accordingly,  the  air  is  introduced  parobol ical ly  into  the  plume 
as  shown  in  Figure  5.  The  crossover  point  defines  an  envelope 
within  which  there  is  no  admixture  of  the  atmosphere.  At  a  given 
point  above  the  candle  surface,  the  parabolic  equations  can  be 
solved  to  give  the  necessary  boundaries  to  compute  the 
composition  of  the  flame  at  each  point  across  the  plume.  Above 
the  crossover  point,  the  per  cent  of  the  atmosphere  at  a  given 
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depth  in  the  plume  is  given  by  the  sum  of  the  mixing  parameters. 

The  physical  boundaries  for  the  flare  plume  were  established 
as  the  limit  for  mixing  with  90%  air.  There  is  another  boundary 
that  must  be  established  for  the  model  also.  This  is  the  optical 
boundary.  Since  the  model  is  concerned  with  the  resonance-line 
continuum  of  an  alkali  metal,  the  optical  boundary  will  be 
defined  as  a  point  in  the  plume  where  the  metal  atom 
concentration  has  decreased  to  1%  of  its  value  in  the  center  of 
the  plume.  The  change  in  concentration  across  the  radius  of  the 
plume  is  due  to  dilution  with  the  atmosphere,  changes  in 
temperature  across  the  plume,  and  reaction  of  the  metal  atoms 
with  other  plume  species.  It  may  be  seen  that  the  optical 
boundary  is  very  close  to  the  physical  boundary  but  the 
distinction  is  important,  especially  with  the  lithium  flames. 

The  temperature  at  the  physical  boundary  is  generally  around  1200 
°K  which  is  close  to  the  value  chosen  by  Douda(l).  The 
temperature  at  the  optical  boundary  is  around  1R00  °K  fer  all  the 
metals  except  lithium  where  it  is  around  2000  °K.  These 
temperatures  are  considered  reasonable  for  both  disappearance  of 
the  metal  atoms  and  the  cessation  of  emission. 

Heat  Loss  Correction 

One  other  factor  must  be  accounted  for  to  allow  the  flame 
model  to  simulate  processes  occurring  in  the  flame.  This  factor 
is  heat  loss  along  the  flare  axis.  There  are  losses  of  the  heat 
of  combustion  due  to  convective,  conductive  and  radiative 
processes  occurring  in  the  plume.  Hamrick,  Tanner,  and  Douda 
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(1  ,  7.9, 30,  34 ,  3  5)  have*  made  calculations  and  measurements  along  the 


plume  axis  to  estimate  the  temperature  along  the  axis.  Using 
their  values  and  combining  this  with  the  assumption  that  the 
visible  radiation  limit  is  1800°K  at  d/I',  a  correction  profile 
for  temperature  can  be  calculated.  A  parabola  is  fit  with  the 
computed  adiabatic  temperature  as  the  vertex  on  the  y-axis  and 
the  height  above  the  candle  surface  along  the  x-axis.  The  flame 
temperature  is  set  at  1800°K  at  d/F  on  the  x-axis.  The  flame 
temperature  along  the  centerline  can  then  bo  computed  as  a 
function  of  distance  above  the  candle  surface.  The  maximum 
temperature  in  any  cross-section  is  determined  by  thermochemical 
calculation  with  the  admixture  of  ambient  atmosphere.  This  value 
includes  excess  heat  of  combustion  that  is  stored  in  the 
magnesium  oxide  and  this  tends  to  hold  the  maximum  temperature  to 
less  temperature  drop  than  along  the?  axis.  It  is  assumed  that 
the  maximum  temperature  will  drop  parabol ical ly  by  a  factor  of 

h 

555(d)  K  along  the  length  of  the  visible  plume  axis.  The 

y3 

factor  includes  a  term  for  d  since  mixing  with  the  atmosphere 
in  the  plume  is  not  complete  as  the  diameter  increases.  Using 
these  fixed  parameters  for  the  rate  of  temperature  drop  along  the 
axis,  a  correction  factor  is  calculated  for  each  incremental  step 
in  the  plume  cross  section  to  adjust  the  temperature  for  heat 
loss.  All  modes  of  heat  loss  are  included  in  the  correction 
factor.  The  values  chosen  are  related  to  the  observed  and 
calculated  temperatures  of  other  researchers  and  will  be  shown  in 
later  sections  to  give  good  agreement  between  calculated  and 
measured  spectral  profiles.  The  actual  equations  used  for  the 
temperature  correction  are  given  in  Appendix  B. 


THE  COMPUTER  PROGRAM 


All  of  the  parameters  required  to  solve  the  radiative 
transfer  equation  may  now  be  computed  from  the  mathematical  model 
described.  To  enable  the  solution  to  be  run  on  a  microcomputer, 
several  tables  of  values  were  calculated  and  stored  as  disk 
files.  When  the  desired  metal  and  formulation  have  been 
selected,  the  appropriate  tables  of  values  are  loaded  into 
memory.  Complete  listings  of  the  main  program  and  the  special 
programs  used  to  construct  the  tables  are  given  in  Appendix  C. 

The  Voigt  function  values  were  computed  using  a  set  of 
parameter  values  in  a  Fortran  program  written  by  G.  Rybicki  and 
included  in  Douda's  program.  It  was  found  that  the  log  of  the 
Voigt  function  at  any  given  wavelength  is  linear  with  respect  to 
the  log  of  the  Voigt  function  «-  parameter.  By  using  the  slope 
and  intercept  of  the  log  function  line,  the  Voigt  function  can  be 
readily  calculated  for  any  Voigt  function  cl  parameter  using  the 
microcomputer  program.  The  Voigt  functions  were  taken  for  five 
different  cl  parameter  values  and  the  slope  and  intercept  of  the 
log  values  computed  by  the  method  of  least  squares.  The  tables 
for  use  in  the  new  model  were  constructed  for  each  metal  by 
selecting  representative  wavelengths,  calculating  the  slope  and 
intercept  of  the  log  values  and  storing  the  values  in  a  disk 
file.  The  actual  tables  used  are  shown  in  Appendix  C  as  Tables 
C-l  through  C-5.  Plots  of  the  slope  and  intercept  as  a  function 
of  wavelength  are  shown  in  Figures  C-l  through  C-5.  The 
discontinuity  of  the  functions  at  the  resonance  wavelengths  is 


apparent 
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To  enable  the  calculation  of  temperature  and  specie 
concentration  through  a  cross  section  of  the  plume,  a  depth 
index,  Z,  is  used.  The  plume  is  divided  into  0.1  cm.  segments 
and  the  position  at  the  plume  surface  is  identified  as  Z=l.  At 
each  step  into  the  plume  of  0.1  cm.,  the  Z  value  is  incremented 
so  that  Z  times  0.1  will  give  the  depth  of  the  plume  in  cm.  The 
per  cent  air  that  has  mixed  with  the  combustion  products  is 
calculated  and  the  temperature,  Voigt  function  ci  parameter, 
total  gas  concentration,  and  metal  gas  concentration  are 
calculated  by  interpolation  from  the  appropriate  look-up  table. 
The  gas  concentrations  and  the  Voigt  function  a.  parameter  are 
adjusted  for  the  chosen  pressure.  The  temperature  profile  is 
next  corrected  for  heat  loss  and  the  Voigt  function  a.  parameter 
is  corrected  for  the  temperature  change  at  each  point  in  the 
plume.  Once  the  temperature,  concentration  and  Voigt  a. 
parameter  functions  have  been  calculated  for  each  position  up  to 
the  centerline  of  the  plume,  a  mirror  image  is  constructed  for 
the  back  half  of  the  plume.  The  optical  thickness  is  computed 
for  each  of  the  Z  index  positions  by  a  numerical  integration 
using  the  trapezoidal  rule.  See  Appendix  B.  for  details.  All  of 
the  parameters  are  now  available  for  a  numerical  integration  of 
Eq .  ( 2 )  by  Simpson's  Rule  (2).  The  values  computed  for  each 
wavelength  are  stored  as  a  disk  file  to  be  recalled  and  plotted 
by  a  special  plotting  routine. 
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EXPERIMENTAL  RESULTS 
Computed  Radiation  and  Temperature  Profiles 

The  mathematical  model  of  the  flare  plume  results  in  a 
predicted  double  peak  in  the  temperature  profile  of  fuel  rich 
compositions.  See  Figure  6.  This  is  due  to  the  reaction  of  the 
excess  fuel  with  the  oxygen  in  the  air  to  maximize  the 
temperature  at  a  radial  distance  that  is  a  function  of  the 
formulation.  This  double  peak  has  been  determined  to  be  a  real 
feature  of  flames  similar  to  these  pyrotechnic  flares 
(9,22,37,38).  This  is  considered  additional  justification  for 
the  selection  of  the  mathematical  model  as  representing  the 
average  flame  processes. 

In  the  solution  of  Eq. (2)  for  the  radiant  power  output,  the 
radiant  power  at  each  point  in  the  flare  plume  is  calculated. 

The  curves  shown  in  Figures  7  through  11  show  the  calculated 
radiant  power  for  five  different  wavelengths  of  a  sodium 
spectrum.  In  the  wings  of  the  spectral  profile  where  the 
resonance  is  weak,  the  radiant  power  is  proportional  to  the 
temperature  profile.  Little  of  the  radiation  from  the  back  of 
the  plume  is  absorbed  as  it  traverses  the  plume.  As  the 
wavelength  moves  closer  to  the  resonance  wavelength,  the  radiant 
power  output  becomes  highly  distorted  as  the  radiant  energy  from 
the  deeper  areas  of  the  plume  is  strongly  absorbed  by  the  cooler 
sodium  near  the  front  surface.  Very  near  the  resonance 
wavelengths,  only  the  radiant  energy  from  the  outer  few  tenths  of 
a  centimeter  can  escape  from  the  plume.  This  is  the  reason  that 
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the  early  model';  had  to  adjust  the  spectral  line  profile 
coefficient  to  a  much  lower  value  to  get  agreement  with 
experiments.  The  Voigt  function  a  parameter  is  only  about  1 . A 
near  the  surface  where  all  the  radiant  power  is  derived  in  the 
wavelength  regions  near  resonance.  This  is  in  good  agreement 
with  the  value  of  1.2  used  by  Douda  as  the  best  empirical  fit  to 
his  d  a t  a . 

Computed  Spectra 

The  computer  program  has  generated  spectra  for  each  of  the 
alkali  metals  at  760  Tor r  for  each  of  three  formulations 

f 

representing  near  stoichiometry  and  two  cases  of  increasing  fuel 
rich  compositions.  These  spectra  are  shown  as  Figures  A-l 
through  A-15  in  Appendix  A.  All  of  the  spectra  were  computed  at 
a  point  one  flare  diameter  above  the  candle  surface  wi th  a  3  cm. 
diameter  candle.  The  computed  radiant  power  values  were 
normalized  to  give  a  maximum  value  of  1.0  for  the  spectra  and 
show  the  shape  of  the  curves.  Figures  A-16  through  A-25  are 
experimental  spectra  taken  by  Blunt  (36)  and  reproduced  here  for 
purposes  of  comparison. 

Since  a  large  volume  of  sodium  spectra  are  available  from  a 
wide  range  of  experimental  conditions,  an  extensive  comparison  of 
calculated  sodium  spectra  has  been  made.  Figure  A-26  shows  the 
computed  spectra  for  several  heights  above  the  candle  surface. 
Experimental  spectra  taken  by  Blunt  (15,16)  show  a  similar  trend 
of  decreasing  power  and  broadening  as  the  position  moves  up  the 
plume  axis.  A  typical  set  of  experimental  spectra  is  shown  in 
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Figure  A-27  for  comparison.  This  variation  is  due  mainly  to  the 
drop  in  plume  temperature  and  the  change  in  cross-section 
concentrations  due  to  additional  mixing  with  air.  The 
differences  noted  in  these  spectra  point  out  how  essential  the 
definition  of  spatial  resolution  is  when  taking  spectra  of 
pyrotechnic  flames. 

Figures  A-28  through  A-41  show  the  computed  spectra  at  different 

ambient  pressures.  These  spectra  are  compared  with  the 

experimental  spectra  of  Pouda  and  Blunt  taken  at  different 

ambient  pressures  (3,11,17).  The  agreement  of  the  computed 

/ 

spectra  with  the  experimental  spectra  demonstrates  that  the  model 
is  predicting  observed  spectral  changes  with  reasonable  accuracy. 

Figure  A-42  shows  a  comparison  of  the  radiant  power  output 
of  two  formulations  integrated  over  the  full  length  of  their 
plumes  at  each  wavelength.  These  integrals  compare  the  total 
power  output  without  correction  for  the  response  of  the  human  eye 
to  radiation.  The  areas  under  the  curves  for  the  range  of 
wavelengths  calculated  are  0.006667  and  0.01616  wa  t  u./steraci  i an . 
These  are  in  the  proper  relationship  for  observed  candlepower 
variations  with  increasing  magnesium  content. 

The  diameter  of  the  candle  has  an  effect  on  the  spectral 
profile  of  the  metal  also.  Figures  A-43  through  A-46  show  the 
computed  spectra  for  four  different  diameter  candles  with  each 
computed  at  two  candle  diameters  above  the  surface.  The  half- 
width  and  the  peak-to-peak  separation  increase  with  diameter. 

The  changes  are  proportional  to  the  square  root  of  the  diameter 
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of  the  candle.  This  is  in  accord  with  experimental  evidence 
reported  by  Douda  (49)  and  in  a  private  communication  from  D.  C. 
A.  Izod,  Royal  Armament  Research  and  Development  Establishment, 
Fort  Halstead,  U.K. 

The  computation  of  the  spectral  profile  of  lithium  has  been 
made  with  comparison  to  the  results  reported  by  Douda  (10)  and  by 
Taylor  and  Farnell  (31).  The  results  are  shown  in  Figures  A-47 
through  A-52.  The  agreement  with  both  the  shape  and  peak  half¬ 
widths  is  very  good  without  making  any  special  adjustment  of 
parameters . 
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For  general  comparison  purposes,  the  spectra  of  potassium, 
rubidium  and  cesium  were  computed  for  the  cases  of  300  Torr,  76 
Torr  and  6  Torr.  The  peak  separation  and  self-absorption  are 
shown  clearly  in  Figures  A-S3  through  A-61. 
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DISCUSSION 
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One  of  the  important  features  of  the  computer  model  is  the 
simpif ication  of  the  computation  of  the  Voigt  function  for  the 
alkali  metals.  The  method  used  by  Rybicki  (1)  gives  a  solution 
for  the  Voigt  function  at  a  large  number  of  wavelengths  for  a 
single  value  of  the  <x  parameter  but  at  the  expense  of  requiring 
a  very  large  computer  and  considerable  computing  time.  This  does 
not  easily  lend  itself  to  the  calculation  of  the  Voigt  function 
for  varying  u  parameters  through  the  plume  cross-section.  The 
microcomputer  solution  uses  a  limited  but  significant  number  of 
wavelengths  to  cover  both  the  broad  spectral  region  and  a  high 
resolution  region.  By  using  a  correlation  of  Voigt  function  with 
the  Voigt  function  parameter,  the  computation  can  be  rapidly 
effected  on  a  microcomputer  over  a  wide  range  of  formulations  and 
conditions  even  though  the  a.  parameter  varies  at  each  point  in 
the  plume. 

The  broad  range  of  cases  that  has  been  examined  strongly 
supports  the  theoretical  model  that  has  been  developed.  Many  of 
the  computed  spectral  phenomena  have  been  observed  experimentally 
without  a  clear  understanding  of  the  mechanisms  producing  the 
results.  It  is  well  established,  for  instance,  that  an  increase 
in  percent  magnesium  above  the  stoichiometric  level  will  produce 
an  increase  in  candlepower  from  an  illuminant.  It  had  been 
noted,  however,  that  the  stoichiometric  formulation  produced  the 
highest  theoretical  temperature  and  spectral  measurements  showed 
that  the  output  per  unit  area  was  highest  for  the  stoichiometric 


formulation.  It  can  bo  neon  from  the  theoretical  model  that  this 
is  the  predicted  mode  for  these  cases.  The  increase  in  magnesium 
causes  an  increase  in  the  total  surface  area  of  the  flare  plume 
so  that,  even  with  a  decrease  in  the  output  per  unit  area,  the 
overall  result  is  a  net  increase  in  candlepowcr.  The  changes  in 
radiant  power  output  and  the  half-width  of  the  metal  resonance 
lines  with  changes  in  ambient  pressure  or  metal  atom 
concentration  are  accounted  for  in  the  model,  also. 

The  model  accounts  for  a  very  wide  range  of  parameters  in 
the  computation  of  the  resonance  broadening  of  all  the  alkali 
metals.  The  unique  variability  of  lithium  atom  concentration 
with  temperature  has  caused  difficulties  with  previous  models. 

The  present  model  adequately  predicts  both  the  half-width  and  the 
spectral  profile  for  lithium  over  a  wide  range.  This  uniqueness 
is  not  handled  as  a  special  case,  however,  but  is  a  part  of  the 
general  solution. 

In  some  cases  the  experimental  spectra  do  not  always  vary  in 
a  systematic  and  consistent  manner.  This  variation  from  the 
p  .dieted  mode  is  not  viewed  as  a  failure  of  the  model  to  predict 
properly  but  is  a  measure  of  the  difficulties  in  obtaining 
reproducible  spectra  from  a  variable  source  such  as  a  pyrotechnic 
flame.  Combustion  irregularities  and  plume  distortion  due  to 
such  things  as  uneven  case  .  ning  can  cause  the  flare  plume  to 
fluctuate  and  the  region  of  interest  will  be  displaced  from  the 
slit  of  the  spectrograph.  The  spectral  profile  is  sensitive  to 
the  height  above  the  candle  surface  and  also  to  radial 
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displacement  from  the  plume  centerline  at  a  fixed  height  above 
the  surface.  Many  of  the  spectra  reported  in  the  literature  were 
taken  under  less  than  ideal  conditions  and  on  a  variety  of 
instruments.  Much  of  the  spectral  variability  is  due  to  the 
experimental  techniques  used.  The  overall  response  of  the  model 
to  known  shifts  in  performance  correlates  well  with  all  of  the 
observed  trends  that  have  been  established  by  years  of  testing. 
The  quantitative  agreement  with  spectra  that  have  been  taken 
recently  with  closely  documented  parameters  is  remarkable  in  that 
no  ad  hoc  assumptions  are  made  across  the  entire  range  of  alkali 
metals . 

f 

In  summary,  it  has  been  shown  that  a  mathematical 
representation  of  a  py rotcchn i c  flame  can  be  used  to  describe  the 
temperature  and  concentration  profiles  of  any  of  the  alkali 
metals  in  the  flame.  The  computed  values  can  then  be  used  in  a 
two-line  radiative  transfer  model  to  compute  the  spectral  radiant 
power  distribution  for  that  metal.  The  model  has  been  developed 
to  be  useful  on  a  microcomputer  by  using  techniques  to  generate 
the  two-line  Voigt  function  from  simplified  slope  and  intercept 
values  with  the  Voigt  function  a.  parameter. 
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Figure  0. 
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Figure  0-5 
Sodiw  at  760  Torr 
507  Kg,  057.  f!3!!03,  57  tinder 
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Figure  A-6 
Sodiw  at  760  Torr 
607  Mg,  357  Na«03,  57  Binder 
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Figure  A- 7 

Potassium  at  760  Torr 
602  lig,  552  Ki<‘J3,  52  tinder 


Figure  A-8 

Potass iun  at  760  Torr 
502  Kg,  657.  13133,  57.  Binder 


Figure  A-9 

Potasskw  at  760  Torr 
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Figure  A-10 
Rubiditw  at.  7<S0  Torr 
357.  I’A,  607  Ro!-:CGt  57  Pindar 
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Figure  A-12 
Kubidiu*  at  760  Torr 
552  Hg,  102  RbN03»  52  terrier 
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Figure?  A- 17 
Lithiuw  at  630  Torr 
61.6%  Mg,  33.^%  L i N 0 3 ,  5%  Binder 
Reproduced  by  Per«iiision  fron  Ref.  (36) 
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Figure  A- 7.7 
RubidiuM  at  630  Torr 
31,47.  Mg ,  63.6%  RbN03 ,  37.  Binder 
Reproduced  by  Remission  fron  Ref.  (36) 


frAOUNT  !N?£NS1?T  nrATTS/STEPA^tAN/M'CR^N) 


tivtitN-.iK  •  *:c*c>.s 


Figure  A-23 
Rub id ion  at,  630  lorr 
66.0%  Mg,  51.0%  RbN03,  5%  Binder 
Reproduced  by  Remission  iron  Ref.  (36) 
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Figure  A- 2^ 

Cesiun  at  630  Terr 
25.9%  Mg,  69.1%  CsNLKl,  5%  Binder 
Reproduced  by  PerMisoion  f  ron  Be  B .  (36) 
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Figure  A-23 
Sodiv.n  at  760  Tori' 
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Figure  A-  29 
Sodiw  at  760  Torr 
E>cper ir.cfital  Spectrun 
Reproduced  by  f'oi inssion  fro«  Rc-f.(l) 
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Figure  fi-30 
Sodiun  at  630  Torr 
Calculated  Spectre* 
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Figure  A-31 
Sodiun  st  630  lorr 
Lxperinental  Speotrun 
Reproduced  by  Remission  fron  Ref.(l) 
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Figure  A-32 

Sodiun  at  300  Torr 

Sodiun  at  300  Torr 

Experimental  Spectrum 

Calculated  Spec t run 

Reproduced  by  Remission  fron  Ref.(l) 
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Figw e  A-3T 
Sodium  3t  130  Terr 
Calculated  Spectrum 
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Figure  fr-35 
Sodium  at  150  Torr 
Experimental  Spectrum 
Reproduced  by  Permission  from  Ref.(l) 
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Figure  A-36 

Figure  A-37 

Sodium  at  75  Torr 

Sodium  at  75  Torr 

Experimental  Spectrum 

Calculated  Spectrum 

Reproduced  by  Permission  from  Ref.(l) 
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Figure  A-V  , 
Lithium  at  760  Terr 
Calculated  Spectrum 


Fig<rre  A--R8 
Lithium  at  760  Torr 
Expor  mental  Spectrum 
Reproduced  by  remission  fro*  Ref. (10) 
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Figure  A-49 
Lithium  at  760  Torr 
Low  Lithium  Concentration 
Calculated  Spectrum 


Figure  A-5D 
Lithium  at  760  Torr 
Low  Lithium  Concentration 
Experinontal  Spectrum 
Reproduced  by  remission  fron  Ref. (10) 
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APPENDIX  B 


Tempo rat  urn  Co  r  rod  i  o n  Ca  1  cti  1  a t  i  o n s 


The  correction  of  temperature  lorn  is  an  important  factor  in 
the  use  of  the  flame  model.  Figures  13-1  and  B-2  show  graphically 
the  assumptions  made  concerning  heat  loss  along  the  axis  and 
radially  across  the  plume.  The  temperature  along  the  centerline 
of  the  flare  plume  axis  is  defined  by  the  equation: 


TCL  - 


0?.OO-  !.)  .2 


(d/i-y- 


L 


+ 


To 


0) 


where  L  is  the  height,  above  the  candle  surface  and  T0  is  the 
adiabatic  temperature  at  0%  air.  The  definitions  of  d  and  F  are 
the  same  as  discussed  in  the  text.  A  general  solution  for  the 
maximum  temperature  is  defined  as: 


T MAX  -  TO  -  65S  (FL/df  D/s 


(2) 


where  TB  is  the  maximum  temperature  calculated  for  the 
stoichiometric  mixture  with  air.  This  temperature  will  always  be 
near  3000  °K,  the  melting  point  of  magnesium  oxide. 

TPRIME  is  the  calculated  temperature  at  the  centerline 
determined  from  the  interpolated  values  in  the  table  of  values 
for  each  Z  index  position.  TJ  is  the  interpolated  maximum 
temperature  at  the  stoichiometric  air  ratio  for  the  given 
composition.  Now,  to  compute  the  temperature  correction  factors, 
the  following  equations  are  used: 
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For  1  <0=  Z  <=  Z MAX , 

Teona  -  TO)  [  l  +  (tmax/tj  -  I)  J  (3) 

where  ZMAX  is  the  index  value  where  the  stoichiometric 
temperature  occurs. 


For  ZMAX  <  Z  <=  (I+l)/2 


f  fj-r  '  V/JAX  1 

f  Tc~  L—  _  Tma\  ” 

{_ 7  t'lM'i f  Tu 

T ‘iuw]  (4) 

TJ  [  ^ 


where  I  is  the  Z  index  value  for  the  back  of  the  plume.  The 
corrected  values  for  T(Z)  are  then  mi r ror-imaged  from  the 
centerline  to  I.  ' 


The  Voigt  function  CL  parameter  values  are  cc  rected  by 
dividing  each  indexed  value  by  the  same  computed  correction 
factors  and  mirror-imaging  the  values  from  the  centerline  to  I 
also . 


Optica1  Thiel:  ness  Calcnlnti  o  n 

The  values  for  the  optical  thickness  ore  determined  by 
calculating  the  optical  density  at  each  index  point  and  storing 
it  in  an  array.  The  optical  density  is 


for  Z  =  2,3,4 . 1  where  7*  -  0  and  HZ  is  defined  as 

H-2  =  l.33Z(,C,XlO~'l*(tAX?  OS)  ((f) 

where  M  is  the  molecular  weight  of  the  metal,  is  the  resonance 
wavelength,  and  OS  is  the  oscillator  strength.  Now,  optical 
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thickness  is 

opc z)  =  op (z-o  +  [orcz)  +  OP(Z-i)]  (l) 

whore  OP(1)=0  and  DC  is  the  depth  increment.  The  values  for 
optical  thickness  are  stored  in  the  computer  program  as  an  array 
for  use  in  the  calculation  of  the  radiant  energy  profile. 
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5  REM : RADTRANS 

10  DIMA  ( 4  0  0 )  ,  PA  { 1 0 , 3 )  ,TT(400)  ,  TAU  (  4  00 , 3)  ,VO(75,2)  ,OP(400)  ,  RV 
(75)  ,  I V  ( 7  5 ) 

: DEFINTJ, K, Z 

15  NA=6.024E23 

sCl=3.71415E24 
s  C2=l . 4  388H8 

:  REM  C 1  =  2  ’•  P 1  *  h  *  c  *  c  *  TO  GIVE  BT  IN  WATTS/CM** 2/MICRON 

16  GOSUB4000 

20  REM  A(400)  TO  THE  ARRAY  FOR  ?,  A  T  R  AT  EACH  Z  INDEX 

21  REM  PA  (9,-1)  IS  THE  ARRAY  FOR  THE  TABLE  CF  VALUES 

22  REM  THE  1st  TERM  REPRESENTS  THE  l  AIR 

23  REM  THE  2nd  TERM  REPRESENTS  TEMP, VOIGT  a, MOLES  GAS, MOLES 
METAL 

24  REM  TAU (400, 3)  IS  AN  ARRAY  HOLDING  THE  INTERPOLATED  VALUE 
S 

25  REM  TT (400)  IS  AN  ARRAY  OF  CALCULATED  OPTICAL  DENSITY  VAL 
UES 

26  REM  OP (400)  IS  THE  OPTICAL  THICKNESS  FOR  EACH  2  INDEX 

27  REM  VO (75, 2 )  IS  AN  ARRAY  FOR  WAVELENGTH , SLOPE  AND  INTLRCE 

PT 

28  REM  THE  VOIGT  FUNCTION  PHI  IS  CALCULATED  FROM  VO (75, 2) 

35  CLS 

: INPUT" ENTER  CANDLE  DIAMETER  IN  CM.";D 
42  PRINT 

SPRINT" FUEL  FRACTION  AT  TEMP.  " ;TF; "  -  " ; F 
SPRINT 
45  DE=.l 

sREM  DE  IS  THE  DELTA  D  FOR  DEPTH  IN  THE  PLUME  (CM) 

50  INPUT" ENTER  HEIGHT  ABOVE  CANDLE  IN  CM.";L 

55  DY=D/2 

s F1=SQR (1+1/F) 
s  PC= (7G0/PR)  ( .  333333 

56  REM  DY  IS  THE  REFERENCE  X  COORDINATE  FOR  THE  CANDLE  EDGE 

60  LM--  (D/F)  *  PC 

s  LX  =  LM/ (FI *SQR (PC) +] )  [2 

61  REM  LM  IS  THE  MAXIMUM  LENGTH  OF  THE  PLUME 

62  REM  LX  IS  THE  CROSS-OVER  HEIGHT  FOR  THE  01  iR  LINES 

63  REM  PR  IS  THE  PRESSURE  IN  TORR 

65  DM-SQR ( PC* (D*L+D*D) ) 

sDl=- (DM/ 2) 

S  D2  =  DY* (F1*SQR (PC) + 1 ) *SQR ( F* L/D/PC) -DY 

66  REM  DM  IS  THE  DIAMETER  OF  THE  PLUME  AT  HEIGHT  L 

67  REM  D1  IS  THE  OUTSIDE  BOUNDARY  FOR  901.  AIR 

68  REM  D2  IS  THE  INSIDE  BOUNDARY  FOR  0%  AIR 

70  I=1NT (DM/DE+. 5) 

71  REM  I  IS  THE  NUMBER  OF  STEPS  THE  PLUME  IS  DIVIDED  INTO 
75  I  FI  =  I  NT  (1/2)  *2TIIENI-I+1 

80  IM=(I+l)/2 


PAGE  81 


: K-INT (D2./DE) 

81  REM  IM  IS  THK  MIDPOINT  INDHX  V  A  IT  IK 

82  RDM  K  IS  THK  INDHX  VALUK  FOR  TI1K  01  AIR  LINE 
8  5  I FL<=LXTHENGOSUB  1 0D0ELSEGOSUB2000 

90  FORZ  =  l'FOI  M 

91  REM  '/MAX  IS  THE  INDEX  VALUE  FOR  THE  MAXIMUM  TEMPERATURE 
95  IFTAU (Z, 0) >TAU (Z-l ,0) THENZMAX=Z 

100  NKXTZ 

112  TPRIME-TAU ( I M , 0) 

:  TJ-TAU  ('/MAX,  0) 

115  TC  L-  (  (TF-  PA  (  0 , 0  )  )  *  ( I./I.M )  f  2  )  +PA  (0,0) 

116  REM  TC L  IS  THE  CENTERLINE  TEMP  CORRECTED  FOR  HEAT  LOSS 

117  REM  TF  IS  THE  TEMPERATURE  LIMIT  FOR  EXCITATION 

120  TMAX-TB- 5 5 5* ( D [ .  333333) * (L/LM)  [ 2 

: Q=INT ( A ( IM) /] 0 ) 

: Ql- (A (IM) -0*10) /10 

:  TPRIME- (PA (0+1 , 0) - TA (0,0) ) *  Ql I  PA  (0,0) 

121  REM  TMAX  IS  THE  MAX  TEMP  CORRECTED  FOR  HEAT  LOSS 

122  REM  TB  IS  THE  MAXIMUM  FLAME  TEMPERATURE 

123  REM  TPRIME  IS  THE , CENTERLINE  TEMP  WITHOUT  HEAT  LOSS  CORI 

ECTION 

125  T 1 -TMAX/TJ 

:  T  2  =  T  C  L/T  PR  I M  E 
130  FORZ  =  2TOZMAX 

:TAU (Z, 0) =TAU (Z, 0) * (1+(T1-1) * ( (Z-l) / (ZMAX-1 ) ) [2) 

132  TAU (Z, 1) =  T  A  U (Z, 1)/ ( 1 + ( T 1 - 1 ) * ( (Z-l )/ (ZMAX-1) ) [2) 

: NEXTZ 

1  35  I FZMAX= I MTHEN 150 

1/0  I FL< LXTHEN I  X- ( IM  +  K) ELSEIX=IM 

1/5  FORZ= ( ZMAX+ 1 ) TO IM 

:TAU  ( Z  ,0)  =TAU  (Z  ,  0)  *  (  (Z-ZMAX)  t  2/  (IM-ZMAX)  [  2*  (T2-T1)  -I-T1 ) 
1/8  TAU ( Z , 1)=TAU(Z, l)/( (Z-ZMAX) (2/(IM-ZMAX) ( 2* (T2-T1 ) +T 1 ) 

:  NEXTZ 

1/9  I  F I  X <  I MT H E N FOR Z  -•  I  XT 0 1 M 
:TAU (Z , 0) =TAU (Z , 0) *T2 
: TAU ( Z , 1 ) =TAU ( Z , 1 ) /T2 
: NEXTZ 

150  FORZ= (IM+l)TOI 
: FORX=  OTO  3 

: TAU ( Z , X) =TAU (I-Z+1,X) 

: NEXTX , Z 

155  117=  1 . 3 3266E-2 Z *MW*  (XLAM  (  3)  *OS 
160  TT ( 0 ) =0 

: FORZ=2TOI 

:TT(Z)=HZ* (TAU (Z , 3) ) 

: NEXTZ 
165  OP(1)=0 

:  FOR7.  =  2TOI 


PAGE  8? 


:0P  ('/.)  : :  0  P  ( Z  - 1 )  H  m;/2‘-  (TT  (Z-l)  HTT  (Z)  ) 

: MRXTZ 

170  GOSUB3000 

171  RE=TAU  (IM,  3)/ 100 

172  FORZ - 2TO 1 0 

: IFTAU  (Z,  3)  >=RET1IENR.7-Z 
:  Z  =  1 1 

: NEXT ELSE NEXT 

173  RPR  I  NT"  CAGE  NO.  ”  ;  C  N?;"  FOR  "  ;  I  ;  "  WITH  PRESSURE  =  "  ;  PR 
TORN  AND  HEIGHT  AEOVE  CANDLE  -  " ; L 

173  FOUR- 1TONN 

:WL-VO (R , 0) 

:SL=VO(K,  1) 

:  I N  =  V 0  ( R  ,  2) 

176  PRINT0512,  "WAVELENGTH  NO.  -  " ;  R 

180  Z  =RJ 

:GOSUP2100 
:  E=OP  (  Z  )  '••'PHI 

181  I F E >  6  OT  H 2  N E  = 6 0 

182  1 F E  < - 6  0  T !  1 K N E  =60  , 

183  TERM  ( 1 )  --BT*  PHI/EXP  (  E) 

185  Z-I-RJ  +  .l 

:  G  OS  U  B  2 1  0  0 
: E=OP ( Z ) *  PHI 

186  I  FE> 60TIIENTERM  ( 2 )  ~  0 

:  GOT 01 9  0 

18  7  TERM (2) =PHI*BT/EXP (OP ( Z ) *  PH  I ) 

ISO  S U M  =  0 

:  FORZ -  (R.7  +  1 )  TO  ( I  -RJ )  STEP2 
: G  OS  U  B  2 i 0  0 
: E=OP ( Z ) *  PH  I 

191  I F  E  >  6  0  T  H  F.  N  E  -  6  0 

192  I  F E < - 6 0 T E; 1 ' N E  =  6 0 

195  SUM=SUM+  PH  I  * BT/EXP  ( E ) 
sNEXTZ 

200  TERM ( 3) =2* SUM 
: SUM=0 

23.0  FORZ-  (R.7  +  2) TO  ( I -R0-1 )  STEP2 
: GOSUB21 00 
: E=OP ( Z ) *  PHI 

211  IFE>60T11ENE  =  60 

212  I F  E  < - 6  0  T  H  E  N  E  =  6  0 

215  SUM=SUM+ PH  I  * BT/EXP ( E ) 

: NEXTZ 

220  TERM (4) ~  4  *  S  UM 

230  IV (R) -DE/3* (TERM ( ] ) +TERM ( 2 ) +TERM ( 3 ) +TERM ( 4 ) ) 

240  PR  I NT0 8 3 2 , VO ( R , 0 )  ,IV(R) 

250  NEXTR 
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260  '0)-0 

>KR  =  1TONN 

1V(R)>IV(0) THEN IV (0)~IV(R) 

265  ■  U 

270  vR-  l'I'ONN 

: RV (R)~TV(R)/IV (0) 

272  LPR1NTVO ( R , 0 ) , IV(R) ,RV(R) 

275  NEXT R 

: CLOSE 3 

280  LINEINPUT" ENTER  THE  NAME  OF  THE  FILE";NM$ 

281  OPEN " 0" , 3 , NM$ 

:  PRINT  it  3  ,  NN 
285  FORQ- 1TONN 

: PRINTS  3, VO (0, 0)  , IV(Q) ,RV (Q) 
jNEXTQ 
: CLOSE 3 

300  RUN" AUTOPLOT 
1000  J-IM+K 
1010  FORZ-  1T0.1 

:A(Z)=(1-.1)*100*( (Z-J) [2/(l-J) [2) 

:  NEXT 

1020  FORX-ZTOIM 
:  A  ( X)  =  0 
:  NEXT 

1030  FORZ~ 1TO ( IM-1 ) 

: A (1M  +  Z ) =A ( IM-Z ) 

:  NEXT 

1035  F0RZ= 1TOIM 
:  FORX-  0TO3 
:0=1NT(A(Z)/10) 

: 01- ( A ( Z ) -0*10) / 1 0 

104  0  TAU(Z,X)=PA(Q,X)  — (PA(Q,X)-PA(0+1  » X)  )  *01 
: NEXTX , Z 

1045  FORZ-1TO (TM-1) 

: FORX=0TO3 

: TAU ( I M  +  Z , X ) ~TAU ( IM-Z  ,  X) 

: NEXTX , Z 
1050  RETURN 
2000  J=IM+K 
2010  FOR'/ ~  1TOJ 

:A(Z)=(1-.1)*100*(  { Z-J )  [ 2/ ( 1-J )  [2) 

:  NEXT 
2020  P=0 

:  FORZ= (2*IM-K)TOIM 
: A ( Z ) =A ( Z ) +A (J-P) 

: P=P+1 
:  NEXT 

2025  FORZ= 1TO ( IM- 1 ) 
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:  A  ( I  M 4  )  --A(IM-Z) 

:  NEXT 

2035  FORZ- 1TOIM 
: FORX--0TO3 
:  Q=  1  K'T  ( A  (  7. )  /  .1  P ) 

:  Ql-  ( A  ( 7. )  -0*  1 0  )  / 1 0 

204  0  TAU  (7  ,  X )  -  PA  ( Q ,  X )  -  (PA  (Q,  X)  -PA  (CH  ]  ,  X)  )  *Q1 
: NEXTX , Z 

2045  FOR  Z r- 1 T O  ( I M  —  1 ) 

: FORX-0TO3 

:  TAU  ( I  H-t  Z  ,  X)  -TAU  ( I  M-Z  ,  X) 

:  NEXT 

20  50  RETURN' 

2099  REM  ROUTINE  TO  CALCULATE  PLANK  FUNCTION  BT 

2100  E  =  E X P  ( C 2 /  (V>  L*TAU  (  Z  ,  0)  )  ) 

2110  BT=C  1  /  (' .'L  [  5*  (E-l)  ) 

2199  REM  ROUTINE  TO  CALCULATE  Pill  FOR  A  GIVEN  VOIGT  a 

2200  PHI- EXP  (SLOPE''' LOG  (TAU  (Z,  1)  )  +IN) 

2220  RET URN 

2999  REM  ROUTINE  TO  LOAD  ARRAY  WITH  WAVELENGTH ,  SLOPE  AND  IN 

TERCEPT  ' 

3000  FI  S~"  PHI  " -*-ID$ 

3005  OPEN" I" , 3, FIS 

: INPUT! 3 , NN 
3015  FORQ- l’i'ONN 
3020  INPUTS,  VO  (0,0) 

: INPUTS  3 , VO (0,1) 

:  INPUT;"  3,  VO  (Q,2) 

: INPUTf  3, CC 
:  NEXT 
3025  CLOSE 

: RETURN 
4000  CLS 

: PRINT" SELECT  THE  METAL  FROM  THE  LIST  BELOW:" 

4005  PRINT 

•.PRINT"  1  -  SODIUM" 

4010  PRINT"  2  -  POTASSIUM" 

4015  PRINT"  3  -  LITHIUM" 

4020  PRINT"  4  -  RUBIDIUM" 

4025  PRINT"  5  -  CESIUM" 

4030  PRINT 

: PRINT" ENTER  THE  NUMBER  OF  YOUR  CHOICE"; 

4035  A $= IN KEYS 

: IFVAL (A$) <10RVAL ( A  $ ) > 5THEN 4 035ELSEI FA$  =  " " THEN  4 0 35 
4040  I D=VAL (AS) 

4042  PRINT 

! PRINT" ENTER  ATMOSPHERIC  PRESSURE  IN  TORR" ; 

: IN PUT PR 


PACE  8  5 


4045  OK  I  DG0T04  100, 4200, 4300, 4400,4500 
4100  I  DO = "  NA" 

:0S=.936 
:OS(l)=.312 
:  OS  ( 2 )  = .  6  2  4 
: MW= 22.9090 
4105  XL ( 1 ) =  589  5 . 92 
: XL ( 2 ) =5889. 95 
: X LAM~  50  9  2 . 9  4 
4110  TF--]  BOO 
4120  G OS tl P, 500  0 
4130  RETURN 
4  20  0  I D  $ - " K " 

: 0S  = 1.023 
:  OS ( 1 )  = . 339 
: OS ( 2 ) = . 60  2 
: MW=  39 .100 
4205  XL (1) -7699.0 
: XL (2) =7  064. 9 
:  XL.A*’  =  7  G8  1.95 
4210  TF=1800  ' 

4215  GOSUB5GOO 
4220  RETURN 
4300  I D$~" LI " 

: OS=. 753 
: OS (1) =. 251 
: OS ( 2) = . 502 
: MW=  6.939 

4305  XL (1) =6707. 912 
:XL(2) =6707.761 
: XLAM=6707 .837 
4310  T F  = 1 9 0 0 
4320  GOSUR  5000 
4330  RETURN 
4400  ID9=" RB" 

: OS= .997 
:  OS  ( 1 )  .  322 
: OS (2) =.675 
: MW=  8  5.40 
4405  X L ( 1 ) =7947.6 
:XL(2) =7000. 3 
: XLAM=7  073.95 
4410  TF= 1 0 0 0 
4415  GOSUP.6000 
4420  RETURN 
4500  I D $  =  " C S " 

: OS= 1.09 
: 05 ( 1 ) = . 362 


i 


32 

:  f-'.W- 1  3  2  .  9  1 
4505  XL(J)-C521.1 
: XL (2)0943. 5 
: XLAM-R7  32 . 3 
4510  T F” 1 R 0 U 
4520  GOSUB7O00 
: RETURN 
5000  CLS 

:  PRINT"  SELECT  A  FORMULATION  FROM  THU  MENU  BELOW:" 
5002  0 $ ~  ]  !'$  t  ’’ NO 3 " 

: PRINT 


5005 

PRINT" 

1 

-  401  MG 

551 

"  ;  J  $ ;  " 

51 

BINDS 

5010 

PRINT" 

2 

-  501  MG 

4  51 

«>  .  1  C*  .  “ 

t  s  t 

51 

n  JNDE 

5015 

PRINT" 

3 

-  601  MG 

3  51 

it  .  -T  (*•  .  11 

/  U  SJ  / 

5 1 

S 1  ’IDE 

5017 

PRINT" 

4 

-  SPECIAL 

CAS  I 

7  It 

5020  PRINT 

:  PRINT"  ENTER  THE  NUMBER  BY  YOUR  CHOICE"; 

5025  CN?=INKEY? 

:  T  KVA!  (CN  ?  )  <  ]  T.  t  C C  )  >  4 T " ” R  5  0 2  5F  T.S  El  FCN$-  "  "THEN 5 02  5 
5030  IFVAI.  (CNC)  =  4THEN  50 50KLSK 5  1  00 

5050  LIN KIN PUT "ENTER  THE  SPECIAL  MIX  IDENTIFIER  FOR  THE  VOIG 
T  FUNCTION" ;NF$ 

5000  GOT 0  5 1 0 5 

510  0  N F ? * " V OICT"  +1  DC+CL’$ 

5105  0 ! ’ E N " I " , 1 ,NF$ 

51)0  FORX- OYOP 


:  FOHY~OTC)3 
;  I N  PUT  !•  1  ,  PA  ( X ,  Y ) 

: NEXTY , X 
5111  TB-  n 

:  "  x = 0 ’  09 

.  .  ■'  PA  (X  ,  0)  >TF'THCNT!5  =  PA  (X,  0) 

:  NEXTXE  LSENEXTX 
5115  FORX-OTOn 

:  V--  PA  (  X  ,  2 )  *  .  C  8  2  0  5  *  PA  ( X  ,  0 )  *  1  C  0  0 

5120  FA (X, 2) - (PR/760) *PA(X, 2) *NA/V 

5121  PA ( X , 3 )  =  ( PR/7  6  0 ) *  PA ( X , 3 ) *NA/V 

5122  NEXTX 

5123  FORX=OTQ9 

: PA (X, 1) -PA(X, 1) * (PR/760) 

5124  NEXTX 

5125  FORX-0TO9 

: I F PA (X, 0) <TFTHENY=X 
:  X=10 

5126  NEXTX 
5130  X-Y-l 

: F~ 1-0 . 10* (Y-(TF-PA(Y,0) )/(PA (X, 0) -PA (Yf0) ) ) 
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5200  CLOSE 

: RETURN 
6000  CEB 

: PRINT" SELECT  A  FORMULATION  FROM  THE  MENU  BELOW: " 


6002 

J0==ID0  +  " 
: PRINT 

..03" 

6  00  5 

PR TNT" 

1  - 

3  5%  MG 

60% 

"  •  i  <;  .  *• 
f  W  V  ; 

.,1  BINDER 

6010 

PRINT" 

o  — 

15?,  MG 

50% 

H  .  7  c-  .  II 

1  l)  v'  1 

5?  BINDER 

6015 

PRINT" 

3  - 

5  5?.  MG 

10  6 

*1  *  1  c-  *  •• 

,  yJ  1 

5%  BINDER 

6017 

6020 

PRINT " 
PRINT 

4  - 

SPECIAL 

cas; 

E" 

: PRINT” ENTER  THE  NUMBER  BY  YOUR  CHOICE"; 

6025  CN  0=1 ’.'KEYS 

:  IFVAL  (CMC)  COEVAL  (CN$)  >  1THEN6  C  2  5EL3EI FCK  $= . THENG025 

6 0  30  I  FVA L  ( C N $  )  =  1 T 1 1 E E ■'  050RLSE61 0 0 

6050  LIN E INPUT"  ENTER  THE  SPECIAL  MIX  IDENTIFIER  FOR  THE  YOIG 
T  FUNCTION" ; NFS 

60  6  0  GOT 06 1 0 5 

610  0  LEO - " VOIGT "  + 1 P?  +  CN $ 

6305  OPEN" I " , 1 ,NF$ 

6110  FOE  X--0  TOO  , 

: FORY-OTG3 

:  INPUT ;U,PA(X,Y) 

: NEXTY , X 

6111  TB~  0 

:  FOR X-'  0TO9 

• IFPA (X, 0) >TBTH  ENTB-PA ( X ,  0) 

:  NEXT  X  E  L.  5  E  N  E  X  T  X 
6115  FORX--OTOO 

:  V=  PA  ( X ,  2 )  *  .  0  F,  2  0  5  *  PA  ( X ,  0 )  *  1  0  0  0 

61  20  PA  f  X  ,  2)  ---  (PR/760)  *PA(X, 2)  *NA/V 

6121  PA ( X , 3 ) “  ( PR/7  6  0 ) *  PA ( X , 3 ) * NA/V 

6122  KEXTX 

6123  FORX-0TO9 

: PA ( X , 1 ) -  PA ( X , 1 ) *  ( PR/7  60) 

6121  NEXTX 

6125  FORX=0TO9 

:  IFF/.  (X,  0)  <TFTHENY  =  X 
:  X=  1 0 

6126  NEXTX 
6130  X ~ Y  — 1 

: F=l-0. 10* (Y- (TF-PA ( Y , 0 ) ) / ( PA ( X , 0 ) -PA ( Y , 0 )  )  ) 

6200  CLOSE 

: RETURN 
7000  CLS 

: PRINT" SELECT  A  FORMULATION  FROM  THE  MENU  BELOW:" 

7002  JS=ir>$  +  "N03" 

: PRINT 


i 


7 


7005 

PRINT" 

1  -  25%  MG 

70%  "  ;  J  $ ; " 

5% 

BINDER" 

7010 

PRINT" 

2-35%  MG 

60%  " ; J 0 ; " 

5% 

BINDER" 

7015 

PRINT" 

3  -  45%  MG 

50%  " ; J 0 ; " 

5% 

BINDER" 

7017 

PRINT" 

4  -  SPECIAL 

CASE" 

7020 

PRINT 

: PRINT" ENTER  THE  NUMBER  BY  YOUR  CHOICE"; 
7025  CN$=INKEY$ 

:  I I'VAL  (CN$ )  <  10RVAL  (CN$)  >  4THEN7  025ELSEI  FCN$: 
70  30  IFVAL  (CNC)  =4THEN7 0 50ELSE7 100 
7050  LINEINPUT" ENTER  THE  SPECIAL  MIX  IDENTIFIER 
T  FUNCTION" ;NF$ 

7060  GOTO6105 

7100  NF$=" VOIGT" +1 D$+CN$ 

7105  OPEN" I " , 1 , NF$ 

7110  FORX=0TO9 

: FORY=0TO3 
: INPUT# 1 ,  PA  ( X,  Y) 

: NEXTY , X 

7111  TB=0 

: FORX=0TO9 

: IFPA (X ,  0) >TBTHENTB=PA (X , 0) 

: NEXTXELSENEXTX 
7115  FORX=0TO9 

:V=PA(X,2) *.08205* PA (X,0)* 1000 

7120  PA ( X , 2 ) ~ (PR/760) *PA (X, 2) *NA/V 

7121  PA  ( X ,  3 )  *:  (PR/760)  *  PA  (X,  3)  *NA/V 

7122  NEXTX 

7123  FORX-0TO9 
:PA(X,1)=PA(X,1)* (PR/760) 

7124  NEXTX 

7125  FORX=0TO9 

: IFPA (X, 0) <TFTHENY=X 
:  X=  1 0 

7126  NEXTX 
7130  X-Y-l 

:F=1-0.10*(Y-(TF-PA(Y,0) )  /  ( PA  (X , 0) -PA (Y , 0) 
7200  CLOSE 

j RETURN 
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" "THEN7025 
FOR  THE  VOIG 


/ 
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10  REM : PROGRAM  TO  CALCULATE  VOIGT  PARAMETERS  FROM  EQUILIBRIUM 
CONSTANTS 
20  CLS 

: DIMVO (9 , 14) ,  N$ ( 1 4 ) ,N ( 14) , SI (1 3) , Q1 ( 1 3) ,M ( 1 3) 

30  PI=3. 1416 

:NA=6. 024E23 
: K=1 . 3804E-16 
: C=2 . 9979  3E18 
40  Kl=3. 1378CE-08 

:K2=7. 15844F.-07 
:K3=5.66584E+06 
50  FORY-~0TO14 

;READN$ (Y) 

:  NEXT 

60  DATA  "TEMP" , "VOIGT  A"," TOTAL  GAS" f "METAL" , "AR" , "CO" , "C02" , 
"H" , " H  2 " , " H20" , "MG" , "N2" , "02" , " OTHER" , "SOLIDS" 

70  F0RY=4T013 

sREADN(Y) , SIGMA (Y) 

:NEXT 

80  DATA  39.94,1E-16,28.01,41E-16,44.01,50E-16,1.008,1E-16,2.0 
16 , 8E-1 6 

90  DATA  18.016,2.2£-16,24.32,1E-16,28.02,21E-16,32.00,34E-16, 
30, 10E-16 
100  CLS 

: IN  PUT" ENTER  THE  CHEMICAL  SYMBOL  OF  THE  METAL" ;CS$ 

110  I FCS $=” NA" THENX= 1ELSEI FCS$  =  " K"THENX=  2ELSEI FCS$  =  " LI " THENX= 

3ELSEIFCS$="RB"THENX=4ELSEIFCS$="CS"THENX=5ELSEGOTO 

100 

120  ONXGOSUB1 000, 1010, 1020,1030, 1040 
130  N ( 3 ) =MW 
140  F0RY=4T01 3 

:Q1(Y)=SIGMA(Y) *SQR (1/N(3)+1/N(Y) )/MU 
: PRINTNS (Y) ,Q1 (Y) 

:  NEXT 

150  FORX= 1T02000 
: NEXTX 
160  CLS 

...  :  FORX=0TO9 
: F0RY=3T01 2 

170  PRINT064 , "THE  VALUES  ARE  TO  BE  ENTERED  FOR  ";X*10;"%  AIR" 
180  POKE16916, 3 

190  PRINT0256, "ENTER  THE  VALUE  FOR  ";N$(Y); 

200  INPUTVO (X,Y) 

:  CLS 
sNEXTY 
210  CLS 

: INPUT" ENTER  THE  VALUE  FOR  SOLIDS" ; VO (X, 14) 

220  CLS 
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'  ; INPUT" ENTER  THE  COMBUSTION  TEMPERATURE" ; VO (X, 0) 

232  SIGMA  (0)  -“-SIGMA  (3)  *VO  (X ,  0 ) /2000 
240  F0RQ=3T01 2 

$VO(X, 2) =VO(X, 2) +VO(X,Q) 
sNEXTQ 
:  CLS 

250  VO{X,13)=l-VO(X,2)-VO(X,14) 

260  VO (X, 2) =  l-VO (X, 14) 

270  CLS 

: INPUT" ENTER  THE  TOTAL  MOLES/KG  OF  ANY  SPECIE" ;TS 
280  PRINT 

: INPUT" ENTER  THE  TOTAL  MOLE  FRACTIONS  OF  THAT  SPECIE" ;TF 
290  TM=TS/TF 

:M=1000/TM 
300  POKE16916 , 0 
:  CLS 

310  REM  ENTER  THE  ROUTINES  TO  CALCULATE  THE  VARIOUS  BROADENIN 

G  PARAMETERS 

320  MU  = (M*MW) / (M+MW) 

330  LN=XLAM[2*A2/ (2*PI*C) 

340  LL= ( X LAMBDA [ 2*SIGMA ( 0 ) /C) * (8*K*NA*VO (X , 0) /PI/MU) [ . 5* (NA* ( 
VO(X,  2)  ~VO(X,  3) ) /VO (X, 2) /82.05/VO (X,0) ) 

350  LR=K1*F*XLAMCDA[3*VO(X,3)/VO(X,2)/VO(X,0) 

360  LD=K2*XLAMBDA* (VO (X, 0) /MW) [.5 
370  SUM=0 

: F0RY=4T01 3 

: SUM=SUM+VO (X, Y) *Q1 (Y) *VO (X , 0) /2000 
: NEXTY 

380  LQ=K3*XLAMRDA[2/(VO(X,0) [.5)/VO(X,  2) *SUM 
390  VO(X, l)=SQR(LOG(2) ) * ( LN+LL+LR+LQ) /LD 
400  PRINT" VOIGT  a  =  ";VO(X,l) 

:FORJ=1T02000 
: NEXT J 

405  LPRINT" %  AIR  FOR  ";CS$;"  =  ";X*10 

: LPR I NT VO ( X , 0 ) ; LN ; LL ; LR ; LQ ; LD ; VO ( X , 1 ) 

410  NEXTX 

420  REM  SAVE  THE  CALCULATED  PARAMETERS  AS  A  FILE  ’VOIGT’ +CHEM 
ICAL  SYMBOL  FOR  METAL  SPECIE 

430  REM  FILE  TO  CONTAIN  THE  TEMPERATURE , VOIGT  a, TOTAL  NO.  OF 
MOLES  AND  MOLES  OF  METAL. 

440  CLS 

: PRINT 
SPRINT 

: PRINT" ENTER  THE  NUMBER  IDENTIFIER  FOR  THE  CASE" 

450  J$=CS$+"N03" 

460  PRINT 

SPRINT"  1  -  40%  MG  55%  5%  BINDER" 

470  PRINT"  2  -  50%  MG  45%  ";J$;"  5%  BINDER" 


I 
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480  PRINT”  3  -  60%  MG  35%  " ;J$;"  5%  BINDER" 

405  PRINT"  4  -  SPECIAL  CASE" 

490  PRINT 

: PRINT" ENTER  THE  NUMBER  BV  THE  CASE  IN  PROGRESS"; 

500  A$=INKEY$ 

:IFVAL(AS) <10RVAL(A$) >4THEN500ELSEIFA$=" "THEN500 
502  IFVAL (A$) =4THENLIMEINPUT" ENTER  THE  SPECIAL  CASE  I  DENT I  PI E 
R" ; ID$ 

:GOTO520 

510  I D$~" VOIGT"+CS$+A$+" ; 1" 

520  OPEN" 0" , 1 , ID$ 

530  FORX=0TO9 

j  FDRY=0TO3 
i PRINTS l,VO{X,Y) 

:NEXTY,X 
540  CLOSE 
600  GOTOIOO 
1000  MW=22 .9098 
:XL=5892.95 
: SIGMA (3)=60E-16 
; MU - 1  t 

: F=. 982 
: A21= . 65E8 
•.RETURN 
1010  MW=39 . 10 

;XL*7681.95 
: SIGMA ( 3) “60E-1 6 
sMU*1.7 
; F=1 . 021 
: A21=. 39E8 
••RETURN 
1020  MW=6 . 939 

;XL*6707.837 
: SIGMA (3)=45E-16 
:MU». 302 
: F=. 753 
; A21=. 37E8 
; RETURN 
1030  MW=0 5.48 

sXL=7873.95 
i SIGMA ( 3) =75E-16 
:MU=3.718 
: F=. 997 
: A2 1= . 36E8 
{RETURN 

1040  MW=132.91 
:XL«8732.5 
: SIGMA (3)=90E-16 
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5  REM  PH  I SLOPE/BAS 
10  CLEAR100 

i DIM  LI (5,1) ,WL{75) ,SL (75) , IN (75) ,CC (75) 

: DEFINTA , N , Q 
:N=1 

: ONERRORGOTOSOOO 
100  CLS 

: PRINT" LINEAR  REGRESSION  -  2  VARIABLES" 

102  INPUT" ENTER  THE  CHEMICAL  SYMBOL  OE  THE  METAL" ;ID$ 

103  LPRINT" VOIGT  FUNCTION  SLOPES  AND  INTERCEPTS  FOR  ";ID$ 

104  LPRINT"  " 

: LPRINT "WAVELENGTH" ; TAB (15) "SLOPE" ; TAB (30) " INTERCEPT" 
B (45) "CORR.  COEFF." 

110  ft" 

120  PRINT 

.•PRINT"  ENTER  -  INDEPENDENT  VARIABLE  FIRST." 

125  F0RX=1T05 

: PRINT" ENTER  POINT  N0.";X;" 

: INPUTLl (X,0) 

126  NEXT 

130  QQ=QQ+1  ,• 

: INPUT" ENTER  THE  WAVELENGTH" ;WL (QQ) 

135  F0RY=1T05 

: PRINT" ENTER  VOIGT  VALUE  NO.  " ; Y ;  " 

: INPUTLl (Y, 1) 

140  NEXT 
145  N=5 

:GOSUB  150 
:GOTO  190 
150  CLS 
155  G=0 

: FOR  A= 1  TO  N 
160  F=A-2*INT( (A-l)/2) 

165  PR I NTT A B (25* (F-l) ) A; ; LI (A, 0) ; " , " ; LI (A, 1) ; 

170  IF  F>=2  THEN  PRINT 

175  IF  A>=30+G*30  THEN  INPUT"PRESS  ENTER  TO  CONTINUE  LIST- 
A$ 

'  :G=G+1 
180  NEXT  A 

: IF  A/2=INT  (A/2)  THEN  PRINT 
185  RETURN 

190  INPUT" DO  YOU  WANT  TO  CHANGE,  DELETE,  OR  ADD  DATA  (Y/N) 
A$ 

195  IF  A$=" Y"  THEN  500 
200  C=2 

:GOTO230 
205  PRINT 

: PRINT" DATA  CAN  BE  ANALYZED  AS  ONE  OF  THE  FOLLOWING:" 
210  PRINTTAB ( 2 5) " 1 .  LINEAR  -  LINEAR" 
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215  PR I NTT AB (25) "2,  LOG  -  LOG" 

220  PRINTTAB (25) " 3 .  LOG  -  LINEAR" 

225  PRINT 

: INPUT" ENTER  THE  NUMBER  BY  THE  DESIRED  CHOICE" ;C 
230  ONCGOTO  250  ,  235  ,  245 
235  F0RQ=1T0N 

: L=LOG ( LI (0,0) ) 

: LI (Q, 0) =L 
:L=LOG (LI (0,1)) 

: LI (0, 1 ) =L 
240  NEXTQ 

:GOTO  250 
245  FORO=1TON 

.  :  L=LOG (LI (Q, 0) ) 

: LI (0, 0) =L 
: NEXTQ 
:GOTO  250 
250  FOR  A= 1  TO  N 

: X1=X1+L1 (A , 0) 

: Y 1=Y 1+L1 ( A , 1 )  , 

j  X2=X2  +  L1 ( A , 0 ) * L 1 (A, 0) 

:Y2=Y2+L1 (A, 1) *L1 (A, 1) 

! Z  =  Z  +  L1 ( A , 0 ) *  LI (A, 1) 

SNEXT  A 
255  X3-X1/N 

; Y3=Y1/N 
sSl=X2-Xl*X3 
; S2=Y2-Y1*Y  3 
: S=Z-X1*Y  3 
260  B=S/S 1 

: B1=Y3-B*X3 
: S3=S2-B*S 

265  IF  N>2  THEN  S4=S3/(N-2) 

270  PRINT"NUMBER  =  ";N 

275  PRINT  "MEAN  OF  X  =  ";X3,"  OF  Y  =  ";Y3 
280  PRINT" SLOPE  =  " ;B,"Y  -  INTERCEPT  =  ";B1 
285  PRINT 

: PRINT" SUM  OF  SQUARES" , "TOTAL" ; Y2 
290  PRINT" MEAN  =  ";Y3*Y1;"  SLOPE  =  ";B*S 

295  PRINT" RESIDUAL  =  ";S3 
300  PRINT 

: PRINT" STANDARD  DEVIATIONS" 

305  PRINT" X=  " ;SQR(S1/(N-1) ) ,"Y  =  " ; SQR (S2/ (N-l) ) 

310  IF  S4<0  THEN  PRINT 
SPRINT 
sGOTO  325 

315  PRINT" ERROR  =  " ; SQR (S4 ) , "Y-BAR  =  ";SQR(S4/N) 

320  PRINT" SLOPE  =  " ;SQR (S4/S1) , "Y-INTER.  =  " ; SQR (S4* ( 1/N+X3* 
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X3/S1 )  ) 

325  PRINT 

:  IF  S4O0  THENPRINT"  F-RATI 0  FOR  ELOPE  =  "  ;B*S/S4 
330  PRINT"CORR.  COEFF.  =  H ; S/SQR (S 1*S2) 

335  X1=0 

:  Y1  =  0 
s  X2=0 
:  Y2=0 
:  Z=0 

338  INPUT"  PRESS  ENTER  TO  CONT.-";A$ 

340  LPRINTWL(QQ)  ;TAB  ( 1 5) D ; TAB ( 30 ) BI ; TAB ( 4 5) S/SQR (S 1*S2) 
342  SL (00) =B 

: IN (QQ) =B1 

: CC (00) =S/SQR (S1*S2) 

344  GOTO480 
355  CLS 

: PRINT0475, "WORKING" 

:MX ( 1 , 1) =0 
: MX ( 1 , 0 ) =0 
360  MN (1,1) =L1  (1,1) 

:MN(1,0)=L1(1,0)' 

: FOR  A= 1  TO  N 

: I F  LI (A, 0) <MN(1,0)  THEN  MN ( 1 , 0) =L1 (A , 0) 

365  IF  LI ( A , I ) <MN  (1,1)  THEN  MN ( 1 , 1 ) =L1 (A , 1 ) 

370  NEXTA 

375  FOR  A= 1  TO  N 

: IF  LI (A , 0) >MX (1,0)  THEN  MX (1,0) ~L1 (A, 0) 

380  IF  LI (A, 1) >MX(1, 1)  THEN  MX ( 1 , 1 ) =L1 ( A, I ) 

385  NEXTA 

390  HI  =  40/ (MN (1,1) -MX  ( 1 r 1 ) ) 

:J1=-MX(1,I) *H1 
395  H2=100/(MX(1,0)-MN(1,0) ) 

: J2=27-MN ( 1 , 0) *H2 
400  CLS 

: FOR  A=0  TO  40 
: SET ( 26 , A) 

! NEXTA 

: FOR  A=0  TO  30  STEP  10 
: SET ( 27 , A) 

: NEXT  A 

405  FOR  A=26  TO  127 
: SET  (A, 41) 

: NEXTA 

: FOR  A=51  TO  126  STEP  25 
:  S  ET  ( A ,  4  0 ) 

: NEXT  A 

410  IFC=1TIIENPRINT0448,"Y"; 

: PRINT093O , " X" ; 
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415  IFC=?.THENPRINT0448  ,  "  LOG  V”; 

: PRINTP930 , " LOG  X"  ; 

420  I  FC=  3TI1KN  PRINTP4  48  ,  "  Y"  ; 

5  PR I NT 9  9  3  0  ,  " LOG  X"; 

425  FOR  A- 1  TO  N 

: SET (LI ( A ,  0 ) *H2+J2,L1 (A, 1) 

: NEXT A 

4  30  PRINT08  37 , " ” ; 

: PRINT  USING  P$;MN(1,1); 

: PRINT? 5 , " " ; 

: PRINT  USING  P$;MX(1,1); 

435  PRINTU908 , " " ; 

SPRINT  USING  P$;HN(1,0); 
s  PRINTP954 , " " ; 

SPRINT  USING  P$ ; MX  ( 1 , 0) ; 

440  FOR  A~27  TO  125 

s  E=  (  ( (A-J2)/H2)*B+B1)*H1+J1 
445  IF  E< 0  THEN  455 
450  SET (A , E) 

455  NEXT  A  , 

4  60  PRINTU  960 SLOP E=" ; B INTER.  =  "; B1 ; "CORR.  COEFF." ;S/SQR(S 
1*S2) 

465  INPUT"  PRESS  ENTER  TO  CONT.-";A$ 

470  K=0 

:W=31+G*30 
sG=G+l 
: RETURN 
475  GOSUB  605 

sFORQ=lTON 
s  L=  L 1 (0,0) 
sLl(Q,0)=Ll (Q,l) 

:L1  (0, 1)-L 
sNEXTO 
jGOSUB  150 
sGOTO  190 

480  ONCGOTO  500  ,  485  ,  495 
485  F0UQ=1T0N 

: L=EXP ( LI (0, 0) ) 
s  LI (Q , 0) =L 
: L=EXP (LI (Q, 1 ) ) 

;L1 (Q, 1} =L 
490  NEXTQ 

sGOTO  500 
495  FORO=1TON 

j  L=EXP ( LI (Q, 0) ) 

:L1(Q,0)=L 
iNEXTQ 
sGOTO  500 


/ 
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500  PRINT" ENTER  A  COMMAND  TO  ADD  -  DELETE  -  CHANGE  -  SWAP  - 
RERUN  -  END" 

505  INPUT" COMMAND" ;B$ 

510  A$=LEFT$ ( B$ , 1 ) 

515  IF  A  $  =  " A "  THEN  550 
520  IF  A$=" D"  THEN  555 
525  IF  A  $  ~ " C "  THEN  580 
530  IF  A  $  =  " S "  THEN  475 
535  IF  A$="R"THENGOSUB  605 
:GOTO  205 

540  IF  A  $  =  "  E  "  THEN  100 

541  IFA$="Q"THEN130 

545  PRINT" ILLEGAL  COMMAND;  ENTER  AGAIN" 

;GOTO  500 
550  N=N+1 

sGOSUR  605 
; GOTO  135 
555  GOSUB  150 

: PRINT" ENTER  PAIR  NUMBER  TO  BE  DELETED"; 

: INPUT  D 
560  GOSUB  605 

565  IF  D>N  THEN  PRINT" PAIR  ";D;"  NOT  FOUND  " 
sGOTO  500 

570  PRINT" PAIR  ";D;M  FOUND  AND  DELETED" 

575  FOR  A=D  TO  N 

:L1 (A,0)=L1(A+1,0) 

:L1 (A, 1)=L1 (A+1,1) 

5  NEXT  A 
:N- N-l 
; GOTO  190 
580  GOSUB  150 

: INPUT" ENTER  THE  DATA  PAIR  NUMBER  TO  BE  CHANGED" ;D 
585  GOSUB  605 

590  IF  D>N  THEN  PRINT" PAIR  " ; D; "  NOT  FOUND-" 

: GOTO  500 

595  PRINT" PAIR  ";D;n  FOUND  -  ENTER  NEW  VALUES  FOR  BOTH"; 

:INPUT  X,Y 
600  LI (D, 0) =X 

:L1 (D, 1) =Y 
sGOTO  190 
605  X1=0 

:  X2=0 
:  X3=0 
:  Y1  =  0 
:  Y2=0 
:  Y3=0 
:  Z=0 
:S=0 


PACK  97 


:  S 1  •-  0 
:  S2=0 
:  S3=0 
:  S4  =  0 
:  B=0 
:  B1  =  0 
610  RETURN 
1000  CLS 

: INPUT" ENTER  THE  NUMBER  OR  WAVELENGTHS" ; NN 
1005  OPEN" O" , 1 , " PHI "+ID0+" : 1" 

:  PRINTi!  1 , NN 
1010  FOROO-1TONN 

:  PRINTj!  1,WL(QQ)  ,SL(QQ)  ,IN  (QQ)  ,CC(QQ) 

.•NEXT 

1020  CLOSE1 
1030  STOP 

2000  OPEN" I " , 1 f " PHI " +ID0+" : 1" 

2010  FORQO- 6T05 5 

:  INPUT (i  1 , WL  (QQ)  ,5L(QQ)  ,  IN  (QQ)  ,CC(00) 

:  NEXT 

2020  CLOSE1  ' 

2030  STOP 

5000  PRINT" ERROR  NO.  ERR/2+1;"  IN  LINE  " ; ERL 
5005  RESUMENEXT 
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